17 #ifndef __deal2__derivative_form_h 18 #define __deal2__derivative_form_h 20 #include <deal.II/base/tensor.h> 40 template <
int order,
int dim,
int spacedim>
137 <<
"Invalid DerivativeForm index " << arg1);
164 template <
int order,
int dim,
int spacedim>
175 template <
int order,
int dim,
int spacedim>
179 Assert( (dim == spacedim) && (order==1),
180 ExcMessage(
"Only allowed for square tensors."));
181 if ((dim == spacedim) && (order==1))
182 for (
unsigned int j=0; j<dim; ++j)
188 template <
int order,
int dim,
int spacedim>
194 for (
unsigned int j=0; j<spacedim; ++j)
201 template <
int order,
int dim,
int spacedim>
206 Assert( (dim == spacedim) && (order==1),
207 ExcMessage(
"Only allowed for square tensors."));
209 if ((dim == spacedim) && (order==1))
210 for (
unsigned int j=0; j<dim; ++j)
218 template <
int order,
int dim,
int spacedim>
223 Assert( (1 == spacedim) && (order==1),
224 ExcMessage(
"Only allowed for spacedim==1 and order==1."));
234 template <
int order,
int dim,
int spacedim>
246 template <
int order,
int dim,
int spacedim>
258 template <
int order,
int dim,
int spacedim>
262 Assert( (1 == spacedim) && (order==1),
271 template <
int order,
int dim,
int spacedim>
275 Assert( (dim == spacedim) && (order==1),
276 ExcMessage(
"Only allowed for square tensors."));
282 if ((dim == spacedim) && (order==1))
283 for (
unsigned int j=0; j<dim; ++j)
292 template <
int order,
int dim,
int spacedim>
301 for (
unsigned int i=0; i<spacedim; ++i)
302 for (
unsigned int j=0; j<dim; ++j)
303 tt[j][i] = (*
this)[i][j];
310 template <
int order,
int dim,
int spacedim>
317 for (
unsigned int i=0; i<spacedim; ++i)
318 for (
unsigned int j=0; j<dim; ++j)
319 dest[i][j] = (*
this)[i] * T[j];
325 template <
int order,
int dim,
int spacedim>
334 return ::determinant(T);
341 for (
unsigned int i=0; i<dim; ++i)
342 for (
unsigned int j=0; j<dim; ++j)
343 G[i][j] = DF_t[i] * DF_t[j];
353 template <
int order,
int dim,
int spacedim>
371 for (
unsigned int i=0; i<dim; ++i)
372 for (
unsigned int j=0; j<dim; ++j)
373 G[i][j] = DF_t[i] * DF_t[j];
382 template <
int order,
int dim,
int spacedim>
406 template <
int spacedim,
int dim>
413 for (
unsigned int i=0; i<spacedim; ++i)
427 template <
int spacedim,
int dim>
435 for (
unsigned int i=0; i<dim; ++i)
447 template <
int spacedim,
int dim>
455 for (
unsigned int i=0; i<spacedim; ++i)
468 template <
int dim,
int spacedim>
479 DEAL_II_NAMESPACE_CLOSE
::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)