Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
field.h
Go to the documentation of this file.
1# ifndef _RHEOLEF_FIELD_H
2# define _RHEOLEF_FIELD_H
3//
4// This file is part of Rheolef.
5//
6// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
7//
8// Rheolef is free software; you can redistribute it and/or modify
9// it under the terms of the GNU General Public License as published by
10// the Free Software Foundation; either version 2 of the License, or
11// (at your option) any later version.
12//
13// Rheolef is distributed in the hope that it will be useful,
14// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16// GNU General Public License for more details.
17//
18// You should have received a copy of the GNU General Public License
19// along with Rheolef; if not, write to the Free Software
20// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
21//
22// =========================================================================
23// AUTHOR: Pierre.Saramito@imag.fr
24// DATE: 23 february 2011
25
26namespace rheolef {
197} // namespace rheolef
198
199#include "rheolef/field_rdof.icc"
200#include "rheolef/field_wdof.icc"
201#include "rheolef/integrate_option.h" // TO_CLEAN with old interface
202
203namespace rheolef {
204
205namespace details {
206
207// forward declarations:
208template <class T, class M> class field_concat_value;
209template <class FieldWdof> class field_wdof_sliced;
210template <class FieldRdof> class field_rdof_sliced_const;
211
212} // namespace details
213
214// ----------------------------------------------------------------------------
215// main class
216// ----------------------------------------------------------------------------
217// [verbatim_field_basic]
218template <class T, class M = rheo_default_memory_model>
219class field_basic: public details::field_wdof_base<field_basic<T,M>> {
220public :
221// typedefs:
222
225 using size_type = std::size_t;
226 using scalar_type = T;
227 using memory_type = M;
229 using geo_type = geo_basic <float_type,memory_type>;
233 using value_type = T; // TODO: run-time dependent, set it to undeterminated<T>
234 using result_type = T; // TODO: run-time dependent, set it to undeterminated<T>
235 class iterator;
236 class const_iterator;
237
238// allocator/deallocator:
239
241
242 explicit field_basic (
243 const space_type& V,
244 const T& init_value = std::numeric_limits<T>::max());
245
246 void resize (
247 const space_type& V,
248 const T& init_value = std::numeric_limits<T>::max());
249
250 // expressions: field_expr or field_lazy are accepted here
251 // TODO: merge with FieldLazy piecewise_poly
252 template <class Expr, class Sfinae =
253 typename std::enable_if<
256 )
258 >::type>
259 field_basic (const Expr& expr);
260
261 field_basic (const std::initializer_list<details::field_concat_value<T,M> >& init_list);
262
263// assignments:
264
266 field_basic<T,M>& operator= (const std::initializer_list<details::field_concat_value<T,M> >& init_list);
267
268// accessors:
269
270 const space_type& get_space() const { return _V; }
271 const geo_type& get_geo() const { return _V.get_geo(); }
272 std::string get_approx() const { return _V.get_approx(); }
273 valued_type valued_tag() const { return _V.valued_tag(); }
274 const std::string& valued() const { return _V.valued(); }
275#ifdef TO_CLEAN
276 std::string name() const { return _V.name(); }
277 bool have_homogeneous_space (space_basic<T,M>& Xh) const { Xh = get_space(); return true; }
278#endif // TO_CLEAN
279
280// accessors & modifiers to unknown & blocked parts:
281
282 const vec<T,M>& u() const { return _u; }
283 const vec<T,M>& b() const { return _b; }
286
287// accessors to extremas:
288
289 T min() const;
290 T max() const;
291 T max_abs() const;
292 T min_abs() const;
293
294// accessors by degrees-of-freedom (dof):
295
296 const distributor& ownership() const { return get_space().ownership(); }
297 const communicator& comm() const { return ownership().comm(); }
298 size_type ndof() const { return ownership().size(); }
299 size_type dis_ndof() const { return ownership().dis_size(); }
300 T& dof (size_type idof);
301 const T& dof (size_type idof) const;
302 const T& dis_dof (size_type dis_idof) const;
303 // write access to non-local dis_idof changes to others procs
305
306 iterator begin_dof();
307 iterator end_dof();
310
311// input/output:
312
314 odiststream& put (odiststream& ops) const;
315 odiststream& put_field (odiststream& ops) const;
316
317// evaluate uh(x) where x is given locally as hat_x in K:
318
319 T dis_evaluate (const point_basic<T>& x, size_type i_comp = 0) const;
320 T operator() (const point_basic<T>& x) const { return dis_evaluate (x,0); }
322
323// [verbatim_field_basic]
324// internals:
326public:
327
328 // evaluate uh(x) where x is given locally as hat_x in K:
329 // requiers to call field::dis_dof_upgrade() before.
330 T evaluate (const geo_element& K, const point_basic<T>& hat_xq, size_type i_comp = 0) const;
331
332 // propagate changed values shared at partition boundaries to others procs
333 template <class SetOp = details::generic_set_op>
334 void dis_dof_update (const SetOp& = SetOp()) const;
335
336 // assignments from expressions
337 template <class Value>
338 typename std::enable_if<
341 >::type
342 operator= (const Value& value)
343 {
344trace_macro("op=(value)...");
345 check_macro (get_space().name() != "", "field=constant : uninitialized field in affectation");
347trace_macro("op=(value) done");
348 return *this;
349 }
350 // TODO: merge with FieldRdof/Lazy && piecewise_polynomial
351 template <class Expr>
352 typename std::enable_if<
354 && ! details::has_field_rdof_interface <Expr>::value
355 && ! details::is_field_expr_v2_constant <Expr>::value
356 && ! details::is_field <Expr>::value
358 >::type
359 operator= (const Expr& expr)
360 {
361trace_macro("op=(expr)...");
362 space_type Xh;
363 check_macro (expr.have_homogeneous_space (Xh),
364 "field = expr; expr should have homogeneous space. HINT: use field = interpolate(Xh, expr)");
365 if (get_space().name() != Xh.name()) resize (Xh, 0);
367trace_macro("op=(expr) done");
368 return *this;
369 }
370 // TODO: FieldRdof && piecewise_polynomial only: filter
371 template<class FieldRdof>
372 typename std::enable_if<
374 && ! details::is_field <FieldRdof>::value
376 >::type
377 operator= (const FieldRdof& rdof)
378 {
379trace_macro("op=(rdof)...");
380 const space_type& Xh = rdof.get_space();
381 if (get_space().name() != Xh.name()) resize (Xh, 0);
383trace_macro("op=(rdof) done");
384 return *this;
385 }
386 // TODO: FieldLazy && piecewise_polynomial only: filter
387 template<class FieldLazy>
388 typename std::enable_if<
392 >::type
393 operator= (const FieldLazy& lazy)
394 {
395trace_macro("op=(lazy)...");
396 const space_type& Xh = lazy.get_space();
397 if (get_space().name() != Xh.name()) resize (Xh, 0);
399trace_macro("op=(lazy) done");
400 return *this;
401 }
402
403#ifdef TO_CLEAN
404 // TODO: move in field_wdof
405 template<class Expr>
406 typename std::enable_if<
409 >::type
410 operator+= (const Expr&);
411
412 // TODO: move in field_wdof
413 template<class Expr>
414 typename std::enable_if<
417 >::type
418 operator-= (const Expr&);
419
420 template<class FieldLazy, class SetPlusOp>
421 typename std::enable_if<
423 ,void
424 >::type
425 convert_from_field_lazy (const FieldLazy& expr, const SetPlusOp& set_plus_op);
426#endif // TO_CLEAN
427
428// accessors by domains and indexes:
429
430 // TODO: why is it necessary to explicit them while its done implicitely in field_wdof_indirect ?
438 { return rdof_base::operator[] (dom_name); }
439
440
441 size_type size() const { return _V.size(); }
450
451 // old interface
452 template <class Expr>
454 const geo_basic<T,M>& dom,
455 const geo_basic<T,M>& band,
456 const band_basic<T,M>& gh,
457 const Expr& expr,
458 const integrate_option& qopt,
459 bool is_on_band);
460 template <class Expr>
461 void do_integrate (
462 const geo_basic<T,M>& dom,
463 const Expr& expr,
464 const integrate_option& iopt);
465 template <class Expr>
466 void do_integrate (
467 const band_basic<T,M>& gh,
468 const Expr& expr,
469 const integrate_option& iopt);
470
471protected:
474
475// data:
477 mutable vec<T,M> _u;
478 mutable vec<T,M> _b;
481// [verbatim_field_basic_cont]
482};
483
484namespace details {
485
486// concepts:
487template<class T, class M>
488struct is_field<field_basic<T,M>> : std::true_type {};
489
490// field class is not an ordinary function/functor : for compose(field,args...) filtering
491template <typename T, typename M> struct function_traits <field_basic<T,M> > {};
492
493template<class T, class M>
495 using size_type = std::size_t;
496 using scalar_type = T;
497 using memory_type = M;
498};
499
500} // namespace details
501
502// [verbatim_field_basic_cont]
503template <class T, class M>
505
506template <class T, class M>
508
509// [verbatim_field]
512// [verbatim_field]
513
515
516// =================================================================================
517// field::iterator
518// =================================================================================
519// TODO: use field_wdof_sliced iterators
520template <class T, class M>
522protected:
523 typedef typename vec<T,M>::iterator data_t; // random acess
524 typedef typename disarray<space_pair_type,M>::const_iterator iter_t; // forward access
525public:
526
527// typedefs:
528
529 typedef std::forward_iterator_tag iterator_category; // TODO: not fully random yet
531 typedef T value_type;
532 typedef T& reference;
533 typedef T* pointer;
534 typedef std::ptrdiff_t difference_type;
535
536// allocators:
537
539 : _blk_dis_iub_iter(),
540 _blk_dis_iub_incr(1),
541 _u(),
542 _b(),
543 _first_iu(0),
544 _first_ib(0)
545 {}
546 iterator (iter_t blk_dis_iub_iter, data_t u, data_t b, size_type first_iu, size_type first_ib)
547 : _blk_dis_iub_iter(blk_dis_iub_iter),
548 _blk_dis_iub_incr(1),
549 _u(u),
550 _b(b),
551 _first_iu(first_iu),
552 _first_ib(first_ib)
553 {}
554 iterator (iter_t blk_dis_iub_iter, size_type blk_dis_iub_incr, data_t u, data_t b, size_type first_iu, size_type first_ib)
555 : _blk_dis_iub_iter(blk_dis_iub_iter),
556 _blk_dis_iub_incr(blk_dis_iub_incr),
557 _u(u),
558 _b(b),
559 _first_iu(first_iu),
560 _first_ib(first_ib)
561 {}
562
563 void set_increment (size_type incr) { _blk_dis_iub_incr = incr; }
564
565// accessors & modifiers:
566
568 bool blk = (*_blk_dis_iub_iter).is_blocked();
569 size_type dis_iub = (*_blk_dis_iub_iter).dis_iub();
570 size_type iub = (!blk) ? dis_iub - _first_iu : dis_iub - _first_ib;
571 return (!blk) ? _u[iub] : _b[iub];
572 }
573 iterator& operator++ () { _blk_dis_iub_iter += _blk_dis_iub_incr; return *this; }
574 iterator operator++ (int) { iterator tmp = *this; operator++(); return tmp; }
575 iterator& operator+= (difference_type n) { _blk_dis_iub_iter += n*_blk_dis_iub_incr; return *this; }
576 iterator operator+ (difference_type n) const { iterator tmp = *this; return tmp += n; }
577 reference operator[] (size_type n) const { return *(*this + n); }
578
579// comparators:
580
581 bool operator== (const iterator& j) const { return _blk_dis_iub_iter == j._blk_dis_iub_iter; }
582 bool operator!= (const iterator& j) const { return ! operator== (j); }
583//protected:
584 friend class const_iterator;
585// data:
591};
592template <class T, class M>
593inline
596{
597 dis_dof_indexes_requires_update();
598 return iterator (_V.data()._idof2blk_dis_iub.begin(), _u.begin(), _b.begin(), _u.ownership().first_index(), _b.ownership().first_index());
599}
600template <class T, class M>
601inline
604{
605 dis_dof_indexes_requires_update();
606 return iterator (_V.data()._idof2blk_dis_iub.end(), _u.begin(), _b.begin(), _u.ownership().first_index(), _b.ownership().first_index());
607}
608// =================================================================================
609// field::const_iterator
610// =================================================================================
611template <class T, class M>
613protected:
616public:
617
618// typedefs:
619
620 typedef std::forward_iterator_tag iterator_category;
622 typedef T value_type;
623 typedef const T& reference;
624 typedef const T* pointer;
625 typedef std::ptrdiff_t difference_type;
626
627// allocators:
628
630 : _blk_dis_iub_iter(),
631 _blk_dis_iub_incr(1),
632 _u(),
633 _b(),
634 _first_iu(0),
635 _first_ib(0)
636 {}
637 const_iterator (iter_t blk_dis_iub_iter, data_t u, data_t b, size_type first_iu, size_type first_ib)
638 : _blk_dis_iub_iter(blk_dis_iub_iter),
639 _blk_dis_iub_incr(1),
640 _u(u), _b(b),
641 _first_iu(first_iu),
642 _first_ib(first_ib)
643 {}
645 : _blk_dis_iub_iter(i._blk_dis_iub_iter),
646 _blk_dis_iub_incr(i._blk_dis_iub_incr),
647 _u(i._u),
648 _b(i._b),
649 _first_iu(i._first_iu),
650 _first_ib(i._first_ib)
651 {}
652
653 void set_increment (size_type incr) { _blk_dis_iub_incr = incr; }
654
655// accessors & modifiers:
656
658 bool blk = (*_blk_dis_iub_iter).is_blocked();
659 size_type dis_iub = (*_blk_dis_iub_iter).dis_iub();
660 size_type iub = (!blk) ? dis_iub - _first_iu : dis_iub - _first_ib;
661 return (!blk) ? _u[iub] : _b[iub];
662 }
663 const_iterator& operator++ () { _blk_dis_iub_iter += _blk_dis_iub_incr; return *this; }
664 const_iterator operator++ (int) { const_iterator tmp = *this; operator++(); return tmp; }
665 const_iterator& operator+= (difference_type n) { _blk_dis_iub_iter += n*_blk_dis_iub_incr; return *this; }
666 const_iterator operator+ (difference_type n) const { const_iterator tmp = *this; return tmp += n; }
667 reference operator[] (size_type n) const { return *(*this + n); }
668
669// comparators:
670
671 bool operator== (const const_iterator& j) const { return _blk_dis_iub_iter == j._blk_dis_iub_iter; }
672 bool operator!= (const const_iterator& j) const { return ! operator== (j); }
673//protected:
674// data:
680};
681template <class T, class M>
682inline
685{
686 dis_dof_indexes_requires_update();
687 return const_iterator (_V.data()._idof2blk_dis_iub.begin(), _u.begin(), _b.begin(), _u.ownership().first_index(), _b.ownership().first_index());
688}
689template <class T, class M>
690inline
693{
694 dis_dof_indexes_requires_update();
695 return const_iterator (_V.data()._idof2blk_dis_iub.end(), _u.begin(), _b.begin(), _u.ownership().first_index(), _b.ownership().first_index());
696}
697// =================================================================================
698// field: inlined
699// =================================================================================
700template <class T, class M>
701inline
703 : _V(),
704 _u(),
705 _b(),
706 _dis_dof_indexes_requires_update(true),
707 _dis_dof_assembly_requires_update(false)
708{
709}
710// TODO: merge Expr with Fieldlazy
711template<class T, class M>
712template<class Expr, class Sfinae>
713inline
715 : _V (),
716 _u (),
717 _b (),
718 _dis_dof_indexes_requires_update(true),
719 _dis_dof_assembly_requires_update(false)
720{
721 operator= (expr);
722}
723template <class T, class M>
724void
726{
727 _dis_dof_indexes_requires_update = true;
728}
729template <class T, class M>
730void
732{
733 _dis_dof_assembly_requires_update = true;
734}
735template <class T, class M>
736inline
737T&
739{
740 dis_dof_indexes_requires_update();
741 bool blk = _V.is_blocked (idof);
742 size_type dis_iub = _V.dis_iub (idof);
743 if (!blk) {
744 size_type iub = dis_iub - _u.ownership().first_index();
745 return _u [iub];
746 } else {
747 size_type iub = dis_iub - _b.ownership().first_index();
748 return _b [iub];
749 }
750}
751template <class T, class M>
752inline
753const T&
755{
756 bool blk = _V.is_blocked (idof);
757 size_type dis_iub = _V.dis_iub (idof);
758 return (!blk) ? _u.dis_at(dis_iub) : _b.dis_at(dis_iub);
759}
760template <class T, class M>
761template <class SetOp>
762void
763field_basic<T,M>::dis_dof_update (const SetOp& set_op) const
764{
765#ifdef _RHEOLEF_HAVE_MPI
766 std::size_t nproc = ownership().comm().size();
767 if (is_distributed<M>::value && nproc > 1) {
768 std::size_t do_assembly = mpi::all_reduce (ownership().comm(), size_t(_dis_dof_assembly_requires_update), std::plus<std::size_t>());
769 std::size_t do_indexes = mpi::all_reduce (ownership().comm(), size_t(_dis_dof_indexes_requires_update), std::plus<std::size_t>());
770 if (do_assembly) {
771 _u.dis_entry_assembly (set_op);
772 _b.dis_entry_assembly (set_op);
773 }
774 if (do_indexes) {
775 _u.set_dis_indexes (_V.ext_iu_set());
776 _b.set_dis_indexes (_V.ext_ib_set());
777 }
778 }
779#endif // _RHEOLEF_HAVE_MPI
780 _dis_dof_indexes_requires_update = false;
781 _dis_dof_assembly_requires_update = false;
782}
783template <class T, class M>
784inline
785T
787{
788 T val = std::numeric_limits<T>::max();
789 for (const_iterator iter = begin_dof(), last = end_dof(); iter != last; iter++) {
790 val = std::min(val, *iter);
791 }
792#ifdef _RHEOLEF_HAVE_MPI
794 val = mpi::all_reduce (comm(), val, mpi::minimum<T>());
795 }
796#endif // _RHEOLEF_HAVE_MPI
797 return val;
798}
799template <class T, class M>
800inline
801T
803{
804 T val = std::numeric_limits<T>::min();
805 for (const_iterator iter = begin_dof(), last = end_dof(); iter != last; iter++) {
806 val = std::max(val, *iter);
807 }
808#ifdef _RHEOLEF_HAVE_MPI
810 val = mpi::all_reduce (comm(), val, mpi::maximum<T>());
811 }
812#endif // _RHEOLEF_HAVE_MPI
813 return val;
814}
815template <class T, class M>
816inline
817T
819{
820 T val = std::numeric_limits<T>::max();
821 for (const_iterator iter = begin_dof(), last = end_dof(); iter != last; iter++) {
822 val = std::min(val, abs(*iter));
823 }
824#ifdef _RHEOLEF_HAVE_MPI
826 val = mpi::all_reduce (comm(), val, mpi::minimum<T>());
827 }
828#endif // _RHEOLEF_HAVE_MPI
829 return val;
830}
831template <class T, class M>
832inline
833T
835{
836 T val = 0;
837 for (const_iterator iter = begin_dof(), last = end_dof(); iter != last; iter++) {
838 val = std::max(val, abs(*iter));
839 }
840#ifdef _RHEOLEF_HAVE_MPI
842 val = mpi::all_reduce (comm(), val, mpi::maximum<T>());
843 }
844#endif // _RHEOLEF_HAVE_MPI
845 return val;
846}
847template <class T, class M>
848inline
851{
852 return uh.get (ips);
853}
854template <class T, class M>
855inline
858{
859 return uh.put (ops);
860}
861// ----------------------------------
862// special interpolate with field arg
863// ----------------------------------
864// re-interpolation of fields and linear field expressions
865// for change of mesh, of approx, ect
866// not truly lazy, but here for the completeness of the interpolate interface
868template<class T, class M>
869inline
870field_basic<T,M>
872{
873 return interpolate (X2h, u1h);
874}
875// interpolate e.g. uh[0] or uh["boundary"], i.e. a FieldRdof
877template <class T, class M, class FieldRdof>
878inline
879typename std::enable_if<
880 details::has_field_rdof_interface<FieldRdof>::value
881 && ! details::is_field<FieldRdof>::value
882 ,field_basic<T,M>
883>::type
884lazy_interpolate (const space_basic<T,M>& Xh, const FieldRdof& uh)
885{
886 return interpolate (Xh, field_basic<T,M>(uh));
887}
888// --------------------------
889// operator=
890// --------------------------
891template <class T, class M>
892inline
893field_basic<T,M>&
895{
896 _V.operator= (x._V);
897 _u.operator= (x._u);
898 _b.operator= (x._b);
899 _dis_dof_indexes_requires_update = x._dis_dof_indexes_requires_update;
900 _dis_dof_assembly_requires_update = x._dis_dof_assembly_requires_update;
901 return *this;
902}
903#ifdef TO_CLEAN
904template<class T, class M>
905template <class Value>
906inline
907typename std::enable_if<
908 details::is_rheolef_arithmetic<Value>::value
910>::type
911field_basic<T,M>::operator= (const Value& value)
912{
913 check_macro (name() != "", "field=constant : uninitialized field in affectation");
914 std::fill (begin_dof(), end_dof(), value);
915 return *this;
916}
917template<class T, class M>
918template<class FieldLazy>
919typename std::enable_if<
920 details::is_field_lazy<FieldLazy>::value
921 ,field_basic<T, M>&
922>::type
923field_basic<T,M>::operator+= (const FieldLazy& expr)
924{
925 const space_basic<T,M>& Xh = expr.get_space();
926 convert_from_field_lazy (expr, details::generic_set_plus_op());
927 return *this;
928}
929template<class T, class M>
930template<class FieldLazy>
931typename std::enable_if<
932 details::is_field_lazy<FieldLazy>::value
933 ,field_basic<T, M>&
934>::type
935field_basic<T,M>::operator-= (const FieldLazy& expr)
936{
937 const space_basic<T,M>& Xh = expr.get_space();
938 convert_from_field_lazy (expr, details::generic_set_minus_op());
939 return *this;
940}
941#endif // TO_CLEAN
942// --------------------------------------
943// io
944// --------------------------------------
945// TODO: one output for all rdof or lazy & one input for all wdof
946template<class FieldLazy>
947typename std::enable_if<
948 details::has_field_lazy_interface<FieldLazy>::value
949 && ! details::is_field<FieldLazy>::value
950 // details::is_field_lazy<FieldLazy>::value
951 ,odiststream&
952>::type
953operator<< (odiststream& out, const FieldLazy& expr)
954{
955 using T = typename FieldLazy::scalar_type;
956 using M = typename FieldLazy::memory_type;
957 field_basic<T,M> uh = expr;
958 return out << uh;
959}
960#ifdef TODO
961// TODO: as: field vh; din >> vh ; uh[i_comp] = vh;
962template<class FieldWdof>
963inline
964idiststream&
965operator>> (odiststream& ids, details::field_wdof_sliced<FieldWdof>& u);
966#endif // TODO
967// --------------------------------------
968// interaction with field_indirect
969// --------------------------------------
970// TODO: use a template FieldRdof and a dof-based convert function
971// => interaction=convert from field_rdof
972#ifdef TO_CLEAN
973template<class T, class M>
974inline
975details::field_wdof_indirect<field_basic<T,M>>
976field_basic<T,M>::operator[] (const geo_type& dom)
977{
978 dis_dof_indexes_requires_update();
979 return details::field_wdof_indirect<field_basic<T,M>> (*this, dom);
980}
981template<class T, class M>
982inline
983details::field_wdof_indirect<field_basic<T,M>>
984field_basic<T,M>::operator[] (std::string dom_name)
985{
986 dis_dof_indexes_requires_update();
987 return operator[] (get_space().get_geo().operator[] (dom_name));
988}
989template<class T, class M>
990inline
991details::field_rdof_indirect_const<field_basic<T,M>>
992field_basic<T,M>::operator[] (const geo_type& dom) const
993{
994 return details::field_rdof_indirect_const<field_basic<T,M>> (*this, dom);
995}
996template<class T, class M>
997inline
998details::field_rdof_indirect_const<field_basic<T,M>>
999field_basic<T,M>::operator[] (std::string dom_name) const
1000{
1001 return operator[] (get_space().get_geo().operator[] (dom_name));
1002}
1003#endif // TO_CLEAN
1004// --------------------------------------
1005// interaction with field_wdof_sliced
1006// --------------------------------------
1007#ifdef TO_CLEAN
1008template<class T, class M>
1009inline
1010details::field_rdof_sliced_const<field_basic<T,M>>
1012{
1013 return details::field_rdof_sliced_const<field_basic<T,M>> (*this, i_comp);
1014}
1015template<class T, class M>
1016inline
1017details::field_wdof_sliced<field_basic<T,M>>
1019{
1020 dis_dof_indexes_requires_update();
1021 return details::field_wdof_sliced<field_basic<T,M>> (*this, i_comp);
1022}
1023// --------------------------------------
1024// convert from field_lazy
1025// --------------------------------------
1026template<class T, class M>
1027template<class FieldLazy, class SetPlusOp>
1028inline
1029typename std::enable_if<
1030 details::is_field_lazy<FieldLazy>::value
1031 ,void
1032>::type
1033field_basic<T,M>::convert_from_field_lazy (const FieldLazy& expr, const SetPlusOp& set_plus_op)
1034{
1035 details::convert_lazy2wdof (expr, set_plus_op, *this);
1036}
1037#endif // TO_CLEAN
1038
1039}// namespace rheolef
1040# endif /* _RHEOLEF_FIELD_H */
field::size_type size_type
Definition branch.cc:430
field gh(Float epsilon, Float t, const field &uh, const test &v)
see the band page for the full documentation
field_rdof_sliced_const< field_basic< T, M > > operator()(size_type i_comp, size_type j_comp) const
field_rdof_indirect_const< field_basic< T, M > > operator[](const std::string &dom_name) const
std::enable_if< details::is_rheolef_arithmetic< Value >::value, field_wdof_base< field_basic< T, M > > & >::type operator=(const Value &)
field_wdof_sliced< field_basic< T, M > > operator()(size_type i_comp, size_type j_comp)
field_wdof_indirect< field_basic< T, M > > operator[](const std::string &dom_name)
rep::base::const_iterator const_iterator
Definition disarray.h:503
see the distributor page for the full documentation
Definition distributor.h:69
size_type dis_size() const
global and local sizes
size_type size(size_type iproc) const
const communicator_type & comm() const
std::forward_iterator_tag iterator_category
Definition field.h:620
const_iterator(iter_t blk_dis_iub_iter, data_t u, data_t b, size_type first_iu, size_type first_ib)
Definition field.h:637
vec< T, M >::size_type size_type
Definition field.h:621
const_iterator operator++(int)
Definition field.h:664
disarray< space_pair_type, M >::const_iterator iter_t
Definition field.h:615
void set_increment(size_type incr)
Definition field.h:653
vec< T, M >::const_iterator data_t
Definition field.h:614
std::forward_iterator_tag iterator_category
Definition field.h:529
vec< T, M >::size_type size_type
Definition field.h:530
iterator(iter_t blk_dis_iub_iter, data_t u, data_t b, size_type first_iu, size_type first_ib)
Definition field.h:546
disarray< space_pair_type, M >::const_iterator iter_t
Definition field.h:524
void set_increment(size_type incr)
Definition field.h:563
std::ptrdiff_t difference_type
Definition field.h:534
vec< T, M >::iterator data_t
Definition field.h:523
iterator(iter_t blk_dis_iub_iter, size_type blk_dis_iub_incr, data_t u, data_t b, size_type first_iu, size_type first_ib)
Definition field.h:554
void dis_dof_indexes_requires_update() const
Definition field.h:725
T max_abs() const
Definition field.h:834
const T & dis_dof(size_type dis_idof) const
Definition field.cc:377
vec< T, M > _u
Definition field.h:477
vec< T, M > _b
Definition field.h:478
const_iterator begin_dof() const
Definition field.h:684
size_type dis_ndof() const
Definition field.h:299
details::field_wdof_indirect< field_basic< T, M > > operator[](const geo_type &dom)
Definition field.h:431
space_type _V
Definition field.h:476
valued_type valued_tag() const
Definition field.h:273
std::string get_approx() const
Definition field.h:272
bool _dis_dof_indexes_requires_update
Definition field.h:479
idiststream & get(idiststream &ips)
Definition field.cc:122
point_basic< T > dis_vector_evaluate(const point_basic< T > &x) const
Definition field.cc:469
size_type size() const
Definition field.h:441
T & dof(size_type idof)
Definition field.h:738
odiststream & put(odiststream &ops) const
Definition field.cc:367
void dis_dof_assembly_requires_update() const
Definition field.h:731
vec< T, M > & set_b()
Definition field.h:285
void dis_dof_update(const SetOp &=SetOp()) const
Definition field.h:763
void do_integrate_internal(const geo_basic< T, M > &dom, const geo_basic< T, M > &band, const band_basic< T, M > &gh, const Expr &expr, const integrate_option &qopt, bool is_on_band)
const vec< T, M > & b() const
Definition field.h:283
field_basic< T, M > & operator=(const field_basic< T, M > &)
Definition field.h:894
dis_reference dis_dof_entry(size_type dis_idof)
Definition field.cc:398
bool _dis_dof_assembly_requires_update
Definition field.h:480
vec< T, M > & set_u()
Definition field.h:284
odiststream & put_field(odiststream &ops) const
Definition field.cc:303
typename vec< scalar_type, memory_type >::dis_reference dis_reference
Definition field.h:231
T operator()(const point_basic< T > &x) const
Definition field.h:320
const T & dof(size_type idof) const
Definition field.h:754
iterator begin_dof()
Definition field.h:595
T min() const
Definition field.h:786
T max() const
Definition field.h:802
T min_abs() const
Definition field.h:818
const geo_type & get_geo() const
Definition field.h:271
int constraint_process_rank() const
Definition field.h:325
const_iterator end_dof() const
Definition field.h:692
const space_type & get_space() const
Definition field.h:270
iterator end_dof()
Definition field.h:603
const vec< T, M > & u() const
Definition field.h:282
size_type ndof() const
Definition field.h:298
const std::string & valued() const
Definition field.h:274
void do_integrate(const geo_basic< T, M > &dom, const Expr &expr, const integrate_option &iopt)
T evaluate(const geo_element &K, const point_basic< T > &hat_xq, size_type i_comp=0) const
Definition field.cc:414
typename float_traits< T >::type float_type
Definition field.h:228
void resize(const space_type &V, const T &init_value=std::numeric_limits< T >::max())
Definition field.cc:49
const distributor & ownership() const
Definition field.h:296
std::size_t size_type
Definition field.h:225
const communicator & comm() const
Definition field.h:297
T dis_evaluate(const point_basic< T > &x, size_type i_comp=0) const
Definition field.cc:439
see the geo_element page for the full documentation
idiststream: see the diststream page for the full documentation
Definition diststream.h:336
see the integrate_option page for the full documentation
odiststream: see the diststream page for the full documentation
Definition diststream.h:137
see the vec page for the full documentation
Definition vec.h:79
base::iterator iterator
Definition vec.h:91
base::size_type size_type
Definition vec.h:86
base::const_iterator const_iterator
Definition vec.h:92
int constraint_process_rank() const
Definition vec_concat.h:91
field_basic< Float, sequential > field_sequential
Definition field.h:514
field_basic< Float > field
see the field page for the full documentation
Definition field.h:511
#define trace_macro(message)
Definition dis_macros.h:111
void get_geo(istream &in, my_geo &omega)
Expr1::float_type T
Definition field_expr.h:230
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
void assign_with_operator(ForwardIterator first, ForwardIterator last, InputIterator iter_rhs, OpAssign op_assign)
std::enable_if< has_field_lazy_interface< FieldLazy >::value &&has_field_wdof_interface< FieldWdof >::value, FieldWdof & >::type convert_lazy2wdof(const FieldLazy &expr0, const SetPlusOp &my_set_plus_op, FieldWdof &uh)
This file is part of Rheolef.
bool operator!=(const heap_allocator< T1 > &lhs, const heap_allocator< T1 > &rhs)
bool operator==(const heap_allocator< T1 > &lhs, const heap_allocator< T1 > &rhs)
std::ostream & operator<<(std::ostream &os, const catchmark &m)
Definition catchmark.h:99
field_basic< T, M > lazy_interpolate(const space_basic< T, M > &X2h, const field_basic< T, M > &u1h)
see the interpolate page for the full documentation
Definition field.h:871
const_iterator begin_dof() const
std::enable_if< details::is_rheolef_arithmetic< U >::value, ad3_basic< T > & >::type operator-=(ad3_basic< T > &a, const U &b)
Definition ad3.h:306
std::enable_if< details::is_rheolef_arithmetic< U >::value, ad3_basic< T > & >::type operator+=(ad3_basic< T > &a, const U &b)
Definition ad3.h:286
bool have_homogeneous_space(space_basic< scalar_type, memory_type > &Vh) const
std::enable_if< details::is_rheolef_arithmetic< U >::value, ad3_basic< T > >::type operator+(const U &a, const ad3_basic< T > &b)
Definition ad3.h:222
std::istream & operator>>(std::istream &is, const catchmark &m)
Definition catchmark.h:88
field_basic< T, M > interpolate(const space_basic< T, M > &V2h, const field_basic< T, M > &u1h)
see the interpolate page for the full documentation
csr< T, sequential > operator*(const T &lambda, const csr< T, sequential > &a)
Definition csr.h:437
Definition leveque.h:25
Expr1::memory_type M