30 #ifndef __OperatorLinearFree_H
31 #define __OperatorLinearFree_H 1
33 #include <feel/feel.hpp>
34 #include <boost/ref.hpp>
35 #include <boost/next_prior.hpp>
36 #include <boost/type_traits.hpp>
37 #include <boost/tuple/tuple.hpp>
38 #if BOOST_VERSION >= 104700
39 #include <boost/math/special_functions/nonfinite_num_facets.hpp>
50 template<
class DomainSpace,
class DualImageSpace,
typename ExprType>
56 typedef ExprType expr_type;
58 typedef typename super::domain_space_type domain_space_type;
59 typedef typename super::dual_image_space_type dual_image_space_type;
60 typedef typename super::domain_space_ptrtype domain_space_ptrtype;
61 typedef typename super::dual_image_space_ptrtype dual_image_space_ptrtype;
64 typedef typename super::backend_ptrtype backend_ptrtype;
67 typedef boost::shared_ptr<matrix_type> matrix_ptrtype;
69 typedef FsFunctionalLinear<DualImageSpace> image_element_type;
70 typedef typename super::domain_element_type domain_element_type;
71 typedef typename super::dual_image_element_type dual_image_element_type;
73 typedef typename super::domain_element_range_type domain_element_range_type;
74 typedef typename super::domain_element_slice_type domain_element_slice_type;
75 typedef typename super::dual_image_element_range_type dual_image_element_range_type;
76 typedef typename super::dual_image_element_slice_type dual_image_element_slice_type;
79 typedef boost::shared_ptr<this_type> this_ptrtype;
87 dual_image_space_ptrtype dualImageSpace,
88 backend_ptrtype backend,
92 super( domainSpace , dualImageSpace , backend ,
false , pattern ),
100 init( domain_space_ptrtype domainSpace,
101 dual_image_space_ptrtype dualImageSpace,
102 backend_ptrtype backend,
103 bool buildMatrix =
false ,
106 this->setDomainSpace( domainSpace );
107 this->setDualImageSpace( dualImageSpace );
119 apply(
const domain_element_type& de,
120 image_element_type& ie )
const
123 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
124 _trial=this->domainSpace(),
125 _pattern=M_pattern );
128 form2( _trial=this->domainSpace(),
129 _test=this->dualImageSpace(),
131 _pattern=M_pattern) = M_expr;
135 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
136 *_v1 = de;_v1->close();
137 vector_ptrtype _v2( M_backend->newVector( _test=ie.space() ) );
138 M_backend->prod( matrix, _v1, _v2 );
139 ie.container() = *_v2;
144 energy(
const typename domain_space_type::element_type & de,
145 const typename dual_image_space_type::element_type & ie )
const
148 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
149 _trial=this->domainSpace(),
150 _pattern=M_pattern );
152 form2( _trial=this->domainSpace(),
153 _test=this->dualImageSpace(),
155 _pattern=M_pattern) = M_expr;
159 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
160 *_v1 = de;_v1->close();
161 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
163 vector_ptrtype _v3( M_backend->newVector( _test=ie.functionSpace() ) );
164 M_backend->prod( matrix, _v1, _v3 );
170 apply(
const typename domain_space_type::element_type & de,
171 typename dual_image_space_type::element_type & ie )
173 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
174 _trial=this->domainSpace(),
175 _pattern=M_pattern );
177 form2( _trial=this->domainSpace(),
178 _test=this->dualImageSpace(),
180 _pattern=M_pattern) = M_expr;
184 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
185 *_v1 = de;_v1->close();
186 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
187 M_backend->prod( matrix, _v1, _v2 );
188 ie.container() = *_v2;
193 apply(
const domain_element_range_type & de,
194 typename dual_image_space_type::element_type & ie )
196 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
197 _trial=this->domainSpace(),
198 _pattern=M_pattern );
200 form2( _trial=this->domainSpace(),
201 _test=this->dualImageSpace(),
203 _pattern=M_pattern) = M_expr;
207 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
208 *_v1 = de;_v1->close();
209 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
210 M_backend->prod( matrix, _v1, _v2 );
211 ie.container() = *_v2;
215 apply(
const typename domain_space_type::element_type & de,
216 dual_image_element_range_type & ie )
219 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
220 _trial=this->domainSpace(),
221 _pattern=M_pattern );
223 form2( _trial=this->domainSpace(),
224 _test=this->dualImageSpace(),
226 _pattern=M_pattern) = M_expr;
230 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
231 *_v1 = de;_v1->close();
232 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
233 M_backend->prod( matrix, _v1, _v2 );
234 ie.container() = *_v2;
238 apply(
const domain_element_range_type & de,
239 dual_image_element_range_type & ie )
242 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
243 _trial=this->domainSpace(),
244 _pattern=M_pattern );
246 form2( _trial=this->domainSpace(),
247 _test=this->dualImageSpace(),
249 _pattern=M_pattern) = M_expr;
253 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
254 *_v1 = de;_v1->close();
255 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
256 M_backend->prod( matrix, _v1, _v2 );
257 ie.container() = *_v2;
261 apply(
const domain_element_slice_type & de,
262 typename dual_image_space_type::element_type & ie )
264 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
265 _trial=this->domainSpace(),
266 _pattern=M_pattern );
268 form2( _trial=this->domainSpace(),
269 _test=this->dualImageSpace(),
271 _pattern=M_pattern) = M_expr;
275 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
276 *_v1 = de;_v1->close();
277 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
278 M_backend->prod( matrix, _v1, _v2 );
279 ie.container() = *_v2;
284 apply(
const typename domain_space_type::element_type & de,
285 dual_image_element_slice_type & ie )
287 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
288 _trial=this->domainSpace(),
289 _pattern=M_pattern );
291 form2( _trial=this->domainSpace(),
292 _test=this->dualImageSpace(),
294 _pattern=M_pattern) = M_expr;
298 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
299 *_v1 = de;_v1->close();
300 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
301 M_backend->prod( matrix, _v1, _v2 );
302 ie.container() = *_v2;
306 apply( domain_element_slice_type de,
307 dual_image_element_slice_type ie )
309 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
310 _trial=this->domainSpace(),
311 _pattern=M_pattern );
313 form2( _trial=this->domainSpace(),
314 _test=this->dualImageSpace(),
316 _pattern=M_pattern) = M_expr;
320 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
321 *_v1 = de;_v1->close();
322 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
323 M_backend->prod( matrix, _v1, _v2 );
324 ie.container() = *_v2;
330 apply(
const domain_element_range_type & de,
331 dual_image_element_slice_type & ie )
333 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
334 _trial=this->domainSpace(),
335 _pattern=M_pattern );
337 form2( _trial=this->domainSpace(),
338 _test=this->dualImageSpace(),
340 _pattern=M_pattern) = M_expr;
344 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
345 *_v1 = de;_v1->close();
346 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
347 M_backend->prod( matrix, _v1, _v2 );
348 ie.container() = *_v2;
352 apply(
const domain_element_slice_type & de,
353 dual_image_element_range_type & ie )
355 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
356 _trial=this->domainSpace(),
357 _pattern=M_pattern );
359 form2( _trial=this->domainSpace(),
360 _test=this->dualImageSpace(),
362 _pattern=M_pattern) = M_expr;
366 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
367 *_v1 = de;_v1->close();
368 vector_ptrtype _v2( M_backend->newVector( _test=ie.functionSpace() ) );
369 M_backend->prod( matrix, _v1, _v2 );
370 ie.container() = *_v2;
376 const image_element_type& ie )
378 auto matrix = M_backend->newMatrix( _test=this->dualImageSpace() ,
379 _trial=this->domainSpace(),
380 _pattern=M_pattern );
382 form2( _trial=this->domainSpace(),
383 _test=this->dualImageSpace(),
385 _pattern=M_pattern) = M_expr;
389 vector_ptrtype _v1( M_backend->newVector( _test=de.functionSpace() ) );
390 vector_ptrtype _v2( M_backend->newVector( _test=ie.space() ) );
391 *_v2 = ie.container();
392 M_backend->solve( matrix, matrix, _v1, _v2 );
398 this_type& operator=( this_type
const& m )
400 M_backend = m.M_backend;
402 M_pattern = m.M_pattern;
408 virtual void matPtr( matrix_ptrtype & matrix )
410 form2( _trial=this->domainSpace(),
411 _test=this->dualImageSpace(),
413 _pattern=M_pattern) = M_expr;
419 backend_ptrtype M_backend;
428 template<
typename Args>
429 struct compute_opLinearFree_return
431 typedef typename boost::remove_reference<typename parameter::binding<Args, tag::domainSpace>::type>::type::element_type domain_space_type;
432 typedef typename boost::remove_reference<typename parameter::binding<Args, tag::imageSpace>::type>::type::element_type image_space_type;
433 typedef typename boost::remove_reference<typename parameter::binding<Args, tag::expr>::type>::type expr_type;
435 typedef OperatorLinearFree<domain_space_type, image_space_type, expr_type> type;
436 typedef boost::shared_ptr<OperatorLinearFree<domain_space_type, image_space_type, expr_type> > ptrtype;
440 BOOST_PARAMETER_FUNCTION(
441 (
typename Feel::detail::compute_opLinearFree_return<Args>::ptrtype ),
445 ( domainSpace, *( boost::is_convertible<mpl::_,boost::shared_ptr<FunctionSpaceBase> > ) )
446 ( imageSpace, *( boost::is_convertible<mpl::_,boost::shared_ptr<FunctionSpaceBase> > ) )
450 ( backend, *, Backend<
typename Feel::detail::compute_opLinearFree_return<Args>::domain_space_type::value_type>::build() )
451 ( pattern, *, (
size_type)Pattern::COUPLED )
456 Feel::detail::ignore_unused_variable_warning( args );
458 typedef typename Feel::detail::compute_opLinearFree_return<Args>::ptrtype oplinfree_ptrtype;
459 return oplinfree_ptrtype (
new oplinfree_type( domainSpace,imageSpace,backend,expr ) );