From: alberto sartori Date: Tue, 10 Nov 2015 12:06:35 +0000 (+0100) Subject: specialized scalar_product symmetric*tensor with sacado X-Git-Tag: v8.4.0-rc2~236^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1402f19d9a3e3e025668cc0480b3941134e4e5e2;p=dealii.git specialized scalar_product symmetric*tensor with sacado --- diff --git a/include/deal.II/base/sacado_product_type.h b/include/deal.II/base/sacado_product_type.h index 5cf6042568..8767d2c257 100644 --- a/include/deal.II/base/sacado_product_type.h +++ b/include/deal.II/base/sacado_product_type.h @@ -72,6 +72,91 @@ struct EnableIfScalar > typedef Sacado::Fad::DFad type; }; + + +/** + * Compute the scalar product $a:b=\sum_{i,j} a_{ij}b_{ij}$ between two + * tensors $a,b$ of rank 2. We don't use operator* for this + * operation since the product between two tensors is usually assumed to be + * the contraction over the last index of the first tensor and the first index + * of the second tensor, for example $(a\cdot b)_{ij}=\sum_k a_{ik}b_{kj}$. + * + * @relates Tensor @relates SymmetricTensor + */ +template +inline +Sacado::Fad::DFad +scalar_product (const SymmetricTensor<2,dim,Sacado::Fad::DFad > &t1, + const Tensor<2,dim,Number> &t2) +{ + Sacado::Fad::DFad s = 0; + for (unsigned int i=0; ioperator* for this + * operation since the product between two tensors is usually assumed to be + * the contraction over the last index of the first tensor and the first index + * of the second tensor, for example $(a\cdot b)_{ij}=\sum_k a_{ik}b_{kj}$. + * + * @relates Tensor @relates SymmetricTensor + */ +template +inline +Sacado::Fad::DFad +scalar_product (const Tensor<2,dim,Number> &t1, + const SymmetricTensor<2,dim,Sacado::Fad::DFad > &t2) +{ + return scalar_product(t2, t1); +} + + +/** + * Compute the scalar product $a:b=\sum_{i,j} a_{ij}b_{ij}$ between two + * tensors $a,b$ of rank 2. We don't use operator* for this + * operation since the product between two tensors is usually assumed to be + * the contraction over the last index of the first tensor and the first index + * of the second tensor, for example $(a\cdot b)_{ij}=\sum_k a_{ik}b_{kj}$. + * + * @relates Tensor @relates SymmetricTensor + */ +template +inline +Sacado::Fad::DFad +scalar_product (const SymmetricTensor<2,dim,Number> &t1, + const Tensor<2,dim,Sacado::Fad::DFad > &t2) +{ + Sacado::Fad::DFad s = 0; + for (unsigned int i=0; ioperator* for this + * operation since the product between two tensors is usually assumed to be + * the contraction over the last index of the first tensor and the first index + * of the second tensor, for example $(a\cdot b)_{ij}=\sum_k a_{ik}b_{kj}$. + * + * @relates Tensor @relates SymmetricTensor + */ +template +inline +Sacado::Fad::DFad +scalar_product (const Tensor<2,dim,Sacado::Fad::DFad > &t1, + const SymmetricTensor<2,dim,Number> &t2) +{ + return scalar_product(t2, t1); +} + DEAL_II_NAMESPACE_CLOSE #endif // DEAL_II_WITH_TRILINOS diff --git a/tests/base/sacado_product_type_06.cc b/tests/base/sacado_product_type_06.cc new file mode 100644 index 0000000000..9e83608de4 --- /dev/null +++ b/tests/base/sacado_product_type_06.cc @@ -0,0 +1,53 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// test scalar_product between tensors and symmetric tensors with Sacado + +#include +#include + +#include +#include + +#include "../tests.h" + + +int main() +{ + typedef Sacado::Fad::DFad Sdouble; + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + + // check product with Tensor<2,dim> + Tensor<2,2,Sdouble> t; + SymmetricTensor<2,2,double> st; + Sdouble a(2,0,7.0); + Sdouble b(2,1,3.0); + + for (unsigned int i=0; i<2; ++i) + for (unsigned int j=0; j<2; ++j) + { + t[i][j] = 2.*a+i*j*b; + st[i][j] = (1.+(i+1)*(j*2)); + } + + deallog << scalar_product(t,st) << std::endl; + deallog << scalar_product(st,t) << std::endl; + +} diff --git a/tests/base/sacado_product_type_06.with_trilinos=on.output b/tests/base/sacado_product_type_06.with_trilinos=on.output new file mode 100644 index 0000000000..84cd33e0dc --- /dev/null +++ b/tests/base/sacado_product_type_06.with_trilinos=on.output @@ -0,0 +1,3 @@ + +DEAL::127.000 [ 16.0000 5.00000 ] +DEAL::127.000 [ 16.0000 5.00000 ] diff --git a/tests/base/sacado_product_type_07.cc b/tests/base/sacado_product_type_07.cc new file mode 100644 index 0000000000..4cb8f944c5 --- /dev/null +++ b/tests/base/sacado_product_type_07.cc @@ -0,0 +1,56 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// test scalar_product between tensors and symmetric tensors with Sacado + +#include +#include + +#include +#include + +#include "../tests.h" + + +int main() +{ + typedef Sacado::Fad::DFad Sdouble; + typedef Sacado::Fad::DFad SSdouble; + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + + // check product with Tensor<2,dim> + Tensor<2,2,SSdouble> t; + SymmetricTensor<2,2,double> st; + SSdouble a(2,0,7.0); + SSdouble b(2,1,3.0); + a.val() = Sdouble(2,0,7.0); + b.val() = Sdouble(2,1,3.0); + + for (unsigned int i=0; i<2; ++i) + for (unsigned int j=0; j<2; ++j) + { + t[i][j] = 2.*a+i*j*b; + st[i][j] = (1.+(i+1)*(j*2)); + } + + deallog << scalar_product(t,st) << std::endl; + deallog << scalar_product(st,t) << std::endl; + +} diff --git a/tests/base/sacado_product_type_07.with_trilinos=on.output b/tests/base/sacado_product_type_07.with_trilinos=on.output new file mode 100644 index 0000000000..4eef31d147 --- /dev/null +++ b/tests/base/sacado_product_type_07.with_trilinos=on.output @@ -0,0 +1,3 @@ + +DEAL::127.000 [ 16.0000 5.00000 ] [ 16.0000 [ ] 5.00000 [ ] ] +DEAL::127.000 [ 16.0000 5.00000 ] [ 16.0000 [ ] 5.00000 [ ] ] diff --git a/tests/base/sacado_product_type_08.cc b/tests/base/sacado_product_type_08.cc new file mode 100644 index 0000000000..5709c1140c --- /dev/null +++ b/tests/base/sacado_product_type_08.cc @@ -0,0 +1,53 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// test scalar_product between tensors and symmetric tensors with Sacado + +#include +#include + +#include +#include + +#include "../tests.h" + + +int main() +{ + typedef Sacado::Fad::DFad Sdouble; + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + + // check product with Tensor<2,dim> + Tensor<2,2,double> t; + SymmetricTensor<2,2,Sdouble> st; + Sdouble a(2,0,7.0); + Sdouble b(2,1,3.0); + + for (unsigned int i=0; i<2; ++i) + for (unsigned int j=0; j<2; ++j) + { + t[i][j] = (1.+(i+1)*(j*2)); + st[i][j] = 2.*a+i*j*b; + } + + deallog << scalar_product(t,st) << std::endl; + deallog << scalar_product(st,t) << std::endl; + +} diff --git a/tests/base/sacado_product_type_08.with_trilinos=on.output b/tests/base/sacado_product_type_08.with_trilinos=on.output new file mode 100644 index 0000000000..7d96063879 --- /dev/null +++ b/tests/base/sacado_product_type_08.with_trilinos=on.output @@ -0,0 +1,3 @@ + +DEAL::155.000 [ 20.0000 5.00000 ] +DEAL::155.000 [ 20.0000 5.00000 ] diff --git a/tests/base/sacado_product_type_09.cc b/tests/base/sacado_product_type_09.cc new file mode 100644 index 0000000000..9335820889 --- /dev/null +++ b/tests/base/sacado_product_type_09.cc @@ -0,0 +1,56 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2015 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE at +// the top level of the deal.II distribution. +// +// --------------------------------------------------------------------- + + +// test scalar_product between tensors and symmetric tensors with Sacado + +#include +#include + +#include +#include + +#include "../tests.h" + + +int main() +{ + typedef Sacado::Fad::DFad Sdouble; + typedef Sacado::Fad::DFad SSdouble; + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.depth_console(0); + deallog.threshold_double(1.e-10); + + + // check product with Tensor<2,dim> + Tensor<2,2,double> t; + SymmetricTensor<2,2,SSdouble> st; + SSdouble a(2,0,7.0); + SSdouble b(2,1,3.0); + a.val() = Sdouble(2,0,7.0); + b.val() = Sdouble(2,1,3.0); + + for (unsigned int i=0; i<2; ++i) + for (unsigned int j=0; j<2; ++j) + { + t[i][j] = (1.+(i+1)*(j*2)); + st[i][j] = 2.*a+i*j*b; + } + + deallog << scalar_product(t,st) << std::endl; + deallog << scalar_product(st,t) << std::endl; + +} diff --git a/tests/base/sacado_product_type_09.with_trilinos=on.output b/tests/base/sacado_product_type_09.with_trilinos=on.output new file mode 100644 index 0000000000..e3d746998d --- /dev/null +++ b/tests/base/sacado_product_type_09.with_trilinos=on.output @@ -0,0 +1,3 @@ + +DEAL::155.000 [ 20.0000 5.00000 ] [ 20.0000 [ ] 5.00000 [ ] ] +DEAL::155.000 [ 20.0000 5.00000 ] [ 20.0000 [ ] 5.00000 [ ] ]