36 #ifndef OPM_JACOBIANSYSTEM_HPP_HEADER
37 #define OPM_JACOBIANSYSTEM_HPP_HEADER
49 namespace ImplicitTransportDefault {
50 template <
class BaseVec>
55 add(
const BaseVec& x, BaseVec& y) {
56 typedef typename BaseVec::value_type VT;
57 ::std::transform(x.begin(), x.end(),
64 template <
class BaseVec>
70 typedef typename BaseVec::value_type VT;
71 ::std::transform(x.begin(), x.end(),
77 template <
class BaseVec>
82 typedef typename BaseVec::value_type VT;
84 ::std::fill(x.begin(), x.end(), VT(0.0));
88 template <
class BaseVec>
93 assign(
const BaseVec& x, BaseVec& y) {
94 ::std::copy(x.begin(), x.end(), y.begin());
98 template <
class Scalar>
100 assign(
const Scalar& a,
const BaseVec& x, BaseVec& y) {
101 ::std::transform(x.begin(), x.end(),
103 [&a](
auto z) { return a * z; });
107 template <
class Matrix>
110 template <
class BaseVec>
113 template <
class Block>
115 assemble(::std::size_t ndof,
120 for (::std::size_t d = 0; d < ndof; ++d) {
121 vec[i*ndof + d] += b[d];
126 template <
class BaseVec>
132 setSize(::std::size_t ndof, ::std::size_t m) {
140 template <
class BaseVec ,
145 enum { Residual = 0, Increment = 1, Solution = 2 };
149 setSize(::std::size_t ndof, ::std::size_t m) {
151 typedef typename ::std::array<BaseVec, 3>::iterator VAI;
153 for (VAI i = vcoll_.begin(), e = vcoll_.end(); i != e; ++i) {
154 VSzSetter<BaseVec>(*i).setSize(ndof, m);
157 for (::std::size_t i = 0; i <
sizeof (vcoll_) /
sizeof (vcoll_[0]); ++i) {
158 VSzSetter<BaseVec>(vcoll_[i]).setSize(ndof, m);
166 VAdd<BaseVec>::add(vcoll_[ Increment ], vcoll_[ Solution ]);
169 template <
class Block>
171 assembleBlock(::std::size_t ndof,
176 assert (ndof == ndof_);
178 VBlkAsm<BaseVec>::assemble(ndof, i, b, vcoll_[ Residual ]);
181 typedef BaseVec vector_type;
183 const vector_type& increment()
const {
return vcoll_[ Increment ]; }
184 const vector_type& residual ()
const {
return vcoll_[ Residual ]; }
185 const vector_type& solution ()
const {
return vcoll_[ Solution ]; }
188 vector_type& writableIncrement() {
return vcoll_[ Increment ]; }
189 vector_type& writableResidual () {
return vcoll_[ Residual ]; }
190 vector_type& writableSolution () {
return vcoll_[ Solution ]; }
193 ::std::size_t ndof_ ;
198 template <
class Matrix>
199 class MatrixBlockAssembler
220 template <
class Matrix ,
221 class NVecCollection>
226 typedef Matrix matrix_type;
227 typedef MatrixBlockAssembler<Matrix> assembler_type;
228 typedef NVecCollection vector_collection_type;
229 typedef typename NVecCollection::vector_type vector_type;
231 assembler_type& matasm() {
return mba_ ; }
232 NVecCollection& vector() {
return sysvec_ ; }
233 const matrix_type& matrix()
const {
return mba_.matrix(); }
234 matrix_type& writableMatrix() {
return mba_.matrix(); }
237 setSize(::std::size_t ndof,
239 ::std::size_t nnz = 0) {
241 mba_ .setSize(ndof, m, m, nnz);
242 sysvec_.setSize(ndof, m );
249 MatrixBlockAssembler<Matrix> mba_ ;
250 NVecCollection sysvec_;
Definition: JacobianSystem.hpp:222
Definition: JacobianSystem.hpp:108
Definition: JacobianSystem.hpp:144
Definition: JacobianSystem.hpp:51
Definition: JacobianSystem.hpp:89
Definition: JacobianSystem.hpp:111
Definition: JacobianSystem.hpp:65
Definition: JacobianSystem.hpp:127
Definition: JacobianSystem.hpp:78
Inverting small matrices.
Definition: ImplicitAssembly.hpp:43