ln_space
Macros | Functions
functions.h File Reference
#include <deal.II/base/exceptions.h>
#include <deal.II/lac/vector.h>
#include <gsl/gsl_integration.h>
#include <gsl/gsl_sf_gamma.h>
#include <cmath>
#include <cstdlib>
#include <numeric>
#include <stdexcept>
#include <string>
#include <utility>
#include <vector>
Include dependency graph for functions.h:
This graph shows which files directly or indirectly include this file:

Go to the source code of this file.

Macros

#define outer_product_sym_H
 

Functions

template<int dim>
Tensor< 4, dim > get_tensor_operator_G (const SymmetricTensor< 2, dim > &Ma, const SymmetricTensor< 2, dim > &Mb)
 
template<int dim>
Tensor< 4, dim > get_tensor_operator_F_right (const SymmetricTensor< 2, dim > &Ma, const SymmetricTensor< 2, dim > &Mb, const SymmetricTensor< 2, dim > &Mc, const SymmetricTensor< 2, dim > &T)
 
template<int dim>
Tensor< 4, dim > get_tensor_operator_F_left (const SymmetricTensor< 2, dim > &Ma, const SymmetricTensor< 2, dim > &Mb, const SymmetricTensor< 2, dim > &Mc, const SymmetricTensor< 2, dim > &T)
 
template<int dim>
SymmetricTensor< 4, dim > outer_product_sym (const SymmetricTensor< 2, dim > &A, const SymmetricTensor< 2, dim > &B)
 
template<int dim>
SymmetricTensor< 4, dim > outer_product_sym (const SymmetricTensor< 2, dim > &A)
 
template<int dim>
SymmetricTensor< 2, dim > outer_product_sym (const Tensor< 1, dim > &A)
 
template<int dim>
bool symmetry_check (Tensor< 2, dim > &tensor)
 
template<int dim>
bool symmetry_check (const Tensor< 4, dim > &temp)
 
template<int dim>
SymmetricTensor< 4, dim > symmetrize (const Tensor< 4, dim > &tensor)
 

Macro Definition Documentation

#define outer_product_sym_H

Function Documentation

template<int dim>
Tensor<4,dim> get_tensor_operator_F_left ( const SymmetricTensor< 2, dim > &  Ma,
const SymmetricTensor< 2, dim > &  Mb,
const SymmetricTensor< 2, dim > &  Mc,
const SymmetricTensor< 2, dim > &  T 
)

Referenced by ln_space< dim >::post_ln().

82  {
83  Tensor<4,dim> tmp;
84 
85  Tensor<2,dim> temp_tensor = contract<1,0>((Tensor<2,dim>)T, (Tensor<2,dim>)Mb);
86  Tensor<2,dim> MaTMb = contract<1,0>((Tensor<2,dim>)Ma,temp_tensor);
87 
88  for(unsigned int i=0; i<dim; ++i)
89  for(unsigned int j=0; j<dim; ++j)
90  for(unsigned int k=0; k<dim; ++k)
91  for(unsigned int l=k; l<dim; ++l)
92  {
93  double tmp_scalar = MaTMb[i][k] * Mc[j][l] + MaTMb[i][l] * Mc[j][k];
94  tmp[i][j][k][l] = tmp_scalar;
95  tmp[i][j][l][k] = tmp_scalar;
96  }
97 
98  return tmp;
99 }
template<int dim>
Tensor<4,dim> get_tensor_operator_F_right ( const SymmetricTensor< 2, dim > &  Ma,
const SymmetricTensor< 2, dim > &  Mb,
const SymmetricTensor< 2, dim > &  Mc,
const SymmetricTensor< 2, dim > &  T 
)

Referenced by ln_space< dim >::post_ln().

57 {
58  Tensor<4,dim> tmp; // has minor symmetry of indices k,l
59 
60  Tensor<2,dim> temp_tensor = contract<1,0>((Tensor<2,dim>)T, (Tensor<2,dim>)Mc);
61  Tensor<2,dim> MbTMc = contract<1,0>((Tensor<2,dim>)Mb,temp_tensor);
62 
63  for(unsigned int i=0; i<dim; ++i)
64  for(unsigned int j=0; j<dim; ++j)
65  for(unsigned int k=0; k<dim; ++k)
66  for(unsigned int l=k; l<dim; ++l)
67  {
68  double tmp_scalar = Ma[i][k] * MbTMc[j][l] + Ma[i][l] * MbTMc[j][k];
69  tmp[i][j][k][l] = tmp_scalar;
70  tmp[i][j][l][k] = tmp_scalar;
71  }
72 
73  return tmp;
74 }
template<int dim>
Tensor<4,dim> get_tensor_operator_G ( const SymmetricTensor< 2, dim > &  Ma,
const SymmetricTensor< 2, dim > &  Mb 
)

Referenced by ln_space< dim >::post_ln().

34 {
35  Tensor<4,dim> tmp; // has minor symmetry of indices k,l
36 
37  for(unsigned int i=0; i<dim; ++i)
38  for(unsigned int j=0; j<dim; ++j)
39  for(unsigned int k=0; k<dim; ++k)
40  for(unsigned int l=k; l<dim; ++l)
41  {
42  double tmp_scalar = Ma[i][k] * Mb[j][l] + Ma[i][l] * Mb[j][k];
43  tmp[i][j][k][l] = tmp_scalar;
44  tmp[i][j][l][k] = tmp_scalar;
45  }
46 
47  return tmp;
48 }
template<int dim>
SymmetricTensor<4,dim> outer_product_sym ( const SymmetricTensor< 2, dim > &  A,
const SymmetricTensor< 2, dim > &  B 
)

Referenced by ln_space< dim >::plastic_right_cauchy_green_AS(), ln_space< dim >::post_ln(), and ln_space< dim >::pre_ln().

111 {
112  SymmetricTensor<4,dim> D;
113  // Special nested for-loop to access only non-symmetric entries of 4th order sym. tensor
114  // ToDo: still not optimal element 1112 and 1211 are both accessed
115  for ( unsigned int i=0; i<dim; ++i )
116  for ( unsigned int j=i; j<dim; ++j )
117  for ( unsigned int k=i; k<dim; ++k )
118  for ( unsigned int l=k; l<dim; ++l )
119  {
120  double tmp = A[i][j] * B[k][l] + B[i][j] * A[k][l];
121  D[i][j][k][l] = tmp;
122  D[k][l][i][j] = tmp;
123  }
124  return D;
125 }
template<int dim>
SymmetricTensor<4,dim> outer_product_sym ( const SymmetricTensor< 2, dim > &  A)
128 {
129  SymmetricTensor<4,dim> D;
130  // Special nested for-loop to access only non-symmetric entries of 4th order sym. tensor
131  // ToDo: still not optimal element 1112 and 1211 are both accessed
132  for ( unsigned int i=0; i<dim; ++i )
133  for ( unsigned int j=i; j<dim; ++j )
134  for ( unsigned int k=i; k<dim; ++k )
135  for ( unsigned int l=k; l<dim; ++l )
136  {
137  double tmp = A[i][j] * A[k][l];
138  D[i][j][k][l] = tmp;
139  D[k][l][i][j] = tmp;
140  }
141  return D;
142 }
template<int dim>
SymmetricTensor<2,dim> outer_product_sym ( const Tensor< 1, dim > &  A)
145 {
146  SymmetricTensor<2,dim> D;
147  // Special nested for-loop to access only non-symmetric entries of 4th order sym. tensor
148  // ToDo: still not optimal element 1112 and 1211 are both accessed
149  for ( unsigned int i=0; i<dim; ++i )
150  for ( unsigned int j=i; j<dim; ++j )
151  D[i][j] = A[i] * A[j];
152 
153  return D;
154 }
template<int dim>
SymmetricTensor<4,dim> symmetrize ( const Tensor< 4, dim > &  tensor)

References symmetry_check().

Referenced by ln_space< dim >::post_ln(), and ln_space< dim >::pre_ln().

209  {
210  SymmetricTensor<4,dim> sym_tensor;
211  Assert(symmetry_check(tensor), ExcMessage("Tensor to symmetrize is not symmetric!"));
212 
213  Tensor<4,dim> temp_tensor;
214  for(unsigned int i=0; i<dim; ++i)
215  for(unsigned int j=0; j<dim; ++j)
216  for(unsigned int k=0; k<dim; ++k)
217  for(unsigned int l=0; l<dim; ++l)
218  temp_tensor[i][j][k][l] = (tensor[i][j][k][l]+tensor[j][i][k][l]+tensor[i][j][l][k]) / 3.0;
219 
220  for(unsigned int i=0; i<dim; ++i)
221  for(unsigned int j=i; j<dim; ++j)
222  for(unsigned int k=0; k<dim; ++k)
223  for(unsigned int l=k; l<dim; ++l)
224  sym_tensor[i][j][k][l] = temp_tensor[i][j][k][l];
225 // sym_tensor[i][j][k][l] = tensor[i][j][k][l];
226 
227  return sym_tensor;
228 }
bool symmetry_check(Tensor< 2, dim > &tensor)
Definition: functions.h:163
template<int dim>
bool symmetry_check ( Tensor< 2, dim > &  tensor)

Referenced by ln_space< dim >::post_ln(), and symmetrize().

164 {
165  for ( unsigned int i=0; i<dim; i++ )
166  for ( unsigned int j=0; j<dim; j++ )
167  if ( i!=j && ( std::abs(tensor[i][j] - tensor[j][i])>1e-12 ) )
168  return false;
169 
170  return true;
171 }
template<int dim>
bool symmetry_check ( const Tensor< 4, dim > &  temp)
173  {
174  for(unsigned int i=0; i<dim; ++i){
175  for(unsigned int j=0; j<dim; ++j){
176  for(unsigned int k=0; k<dim; ++k){
177  for(unsigned int l=0; l<dim; ++l){
178  // absolute difference check
179  if( ( fabs(temp[i][j][k][l] - temp[j][i][k][l]) > 1e-6 )
180  ||
181  ( fabs(temp[i][j][k][l] - temp[i][j][l][k]) > 1e-6 ) )
182  {
183  // relative difference check
184  if( (fabs(fabs(temp[i][j][k][l] - temp[j][i][k][l])/temp[j][i][k][l])> 1e-8)
185  ||
186  (fabs(fabs(temp[i][j][k][l] - temp[i][j][l][k])/temp[i][j][l][k]) >1e-8) )
187  {
188  deallog<< std::endl;
189  deallog<< "Abs not matching: "<<fabs(temp[i][j][k][l] - temp[j][i][k][l])
190  << " | Rel not matching: "<<fabs(temp[i][j][k][l] - temp[j][i][k][l])/temp[j][i][k][l]
191  << " | Abs not matching: "<<fabs(temp[i][j][k][l] - temp[i][j][l][k])
192  << " | Rel not matching: "<<fabs(temp[i][j][k][l] - temp[i][j][l][k])/temp[i][j][l][k]
193  << " | value ijkl: "<< std::scientific << temp[i][j][k][l]
194  << " | value jikl: "<< std::scientific << temp[j][i][k][l]
195  << " | value ijlk: "<< std::scientific << temp[i][j][l][k]
196  << std::endl;
197  return false;
198  }
199  }
200  }
201  }
202  }
203  }
204  return true;
205 }