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 typedef typename BaseVec::value_type VT;
102 ::std::transform(x.begin(), x.end(),
104 ::std::bind2nd(::std::multiplies<VT>(), a));
108 template <
class Matrix>
111 template <
class BaseVec>
114 template <
class Block>
116 assemble(::std::size_t ndof,
121 for (::std::size_t d = 0; d < ndof; ++d) {
122 vec[i*ndof + d] += b[d];
127 template <
class BaseVec>
133 setSize(::std::size_t ndof, ::std::size_t m) {
141 template <
class BaseVec ,
146 enum { Residual = 0, Increment = 1, Solution = 2 };
150 setSize(::std::size_t ndof, ::std::size_t m) {
152 typedef typename ::std::array<BaseVec, 3>::iterator VAI;
154 for (VAI i = vcoll_.begin(), e = vcoll_.end(); i != e; ++i) {
155 VSzSetter<BaseVec>(*i).setSize(ndof, m);
158 for (::std::size_t i = 0; i <
sizeof (vcoll_) /
sizeof (vcoll_[0]); ++i) {
159 VSzSetter<BaseVec>(vcoll_[i]).setSize(ndof, m);
167 VAdd<BaseVec>::add(vcoll_[ Increment ], vcoll_[ Solution ]);
170 template <
class Block>
172 assembleBlock(::std::size_t ndof,
177 assert (ndof == ndof_);
179 VBlkAsm<BaseVec>::assemble(ndof, i, b, vcoll_[ Residual ]);
182 typedef BaseVec vector_type;
184 const vector_type& increment()
const {
return vcoll_[ Increment ]; }
185 const vector_type& residual ()
const {
return vcoll_[ Residual ]; }
186 const vector_type& solution ()
const {
return vcoll_[ Solution ]; }
189 vector_type& writableIncrement() {
return vcoll_[ Increment ]; }
190 vector_type& writableResidual () {
return vcoll_[ Residual ]; }
191 vector_type& writableSolution () {
return vcoll_[ Solution ]; }
194 ::std::size_t ndof_ ;
199 template <
class Matrix>
200 class MatrixBlockAssembler
221 template <
class Matrix ,
222 class NVecCollection>
227 typedef Matrix matrix_type;
228 typedef MatrixBlockAssembler<Matrix> assembler_type;
229 typedef NVecCollection vector_collection_type;
230 typedef typename NVecCollection::vector_type vector_type;
232 assembler_type& matasm() {
return mba_ ; }
233 NVecCollection& vector() {
return sysvec_ ; }
234 const matrix_type& matrix()
const {
return mba_.matrix(); }
235 matrix_type& writableMatrix() {
return mba_.matrix(); }
238 setSize(::std::size_t ndof,
240 ::std::size_t nnz = 0) {
242 mba_ .setSize(ndof, m, m, nnz);
243 sysvec_.setSize(ndof, m );
250 MatrixBlockAssembler<Matrix> mba_ ;
251 NVecCollection sysvec_;
Definition: JacobianSystem.hpp:109
Definition: JacobianSystem.hpp:112
Definition: JacobianSystem.hpp:128
Definition: JacobianSystem.hpp:145
Definition: JacobianSystem.hpp:223
Definition: JacobianSystem.hpp:89
Definition: JacobianSystem.hpp:78
Definition: JacobianSystem.hpp:65
Definition: JacobianSystem.hpp:51