Flow123d  master-f44eb46
field.hh
Go to the documentation of this file.
1 /*!
2  *
3  * Copyright (C) 2015 Technical University of Liberec. All rights reserved.
4  *
5  * This program is free software; you can redistribute it and/or modify it under
6  * the terms of the GNU General Public License version 3 as published by the
7  * Free Software Foundation. (http://www.gnu.org/licenses/gpl-3.0.en.html)
8  *
9  * This program is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12  *
13  *
14  * @file field.hh
15  * @brief
16  */
17 
18 #ifndef FIELD_HH_
19 #define FIELD_HH_
20 
21 #include <stdio.h> // for sprintf
22 #include <string.h> // for memcpy
23 #include <algorithm> // for find, min
24 #include <boost/circular_buffer.hpp>
25 #include <memory> // for dynamic_pointe...
26 #include <new> // for operator new[]
27 #include <ostream> // for basic_ostream:...
28 #include <string> // for basic_string
29 #include <utility> // for pair
30 #include <vector> // for vector
31 #include <armadillo>
32 #include "fields/field_algo_base.hh" // for FieldAlgorithm...
33 #include "fields/field_common.hh" // for FieldCommon::T...
34 #include "fields/field_values.hh" // for FieldValue<>::...
35 #include "fields/field_value_cache.hh" // for FieldValueCache
36 #include "input/accessors.hh" // for ExcTypeMismatch
37 #include "input/accessors_impl.hh" // for Record::opt_val
38 #include "input/factory_impl.hh" // for Factory::create
39 #include "input/input_exception.hh" // for FieldCommon::E...
40 #include "input/storage.hh" // for ExcStorageType...
41 #include "input/type_base.hh" // for Array
42 #include "input/type_generic.hh" // for Instance
43 #include "input/type_record.hh" // for Record::ExcRec...
44 #include "input/input_exception.hh" // for Input::Exception
45 #include "io/output_time.hh" // for OutputTime
46 #include "mesh/elements.h" // for Element::dim
47 #include "mesh/region.hh" // for RegionDB::ExcU...
48 #include "system/asserts.hh" // for Assert, ASSERT
49 #include "system/exceptions.hh" // for ExcAssertMsg::...
50 #include "tools/time_governor.hh" // for TimeStep
51 
52 class Mesh;
53 class Observe;
54 class EvalPoints;
55 class BulkPoint;
56 class SidePoint;
57 class FieldSet;
59 template <int spacedim> class ElementAccessor;
60 template <int spacedim, class Value> class FieldFE;
61 namespace detail
62 {
63  template< typename CALLABLE, typename TUPLE, int INDEX >
64  struct model_cache_item;
65 }
66 
67 using namespace std;
68 namespace IT=Input::Type;
69 
70 /**
71  * @brief Class template representing a field with values dependent on: point, element, and region.
72  *
73  * By "field" we mean a mapping of a a pair (Point, Time) to a @p Value, where
74  * Point is from @p spacedim dimensional ambient space, Time is real number (set by @p set_time method),
75  * and @p Value type representing range of the field, which can be: real scalar, integer scalar (a discrete value),
76  * real vector of fixed (compile time) size, real vector of runtime size, or a matrix of fixed dimensions.
77  * Extensions to vectors or matrices of integers, or to variable tensors are possible. For vector and matrix values
78  * we use classes provided by Armadillo library for linear algebra.
79  * The @p Value template parameter should FieldValue<> template, usual choices are:
80  * FieldValue<spacedim>::Scalar, FieldValue<spacedim>::Integer, FieldValue<spacedim>::Enum,
81  * FieldValue<spacedim>::VectorFixed, FieldValue<spacedim>::TensorFixed.
82  *
83  * This class assign particular fields (instances of descendants of FiledBase) to the regions. It keeps a table of pointers to fields for every possible bulk
84  * region index (very same functionality, but for boundary regions is provided by @p BCField class). This class has interface very similar to FiledBase, however
85  * key methods @p value, and @p value_list are not virtual in this class by contrast these methods are inlined to minimize overhead for
86  * simplest fields like FieldConstant.
87  *
88  * TODO: currently it works only for spacedim==3 since we have only mesh in 3D ambient space.
89  *
90  */
91 template<int spacedim, class Value>
92 class Field : public FieldCommon {
93 public:
94 
96  typedef std::shared_ptr< FieldBaseType > FieldBasePtr;
98  typedef Value ValueType;
99 
100  static const unsigned int space_dim = spacedim;
101 
102 
103  /**
104  * Factory class that creates an instance of FieldBase for field
105  * with name @p field_name based on data in field descriptor @p rec.
106  *
107  * Default implementation in method @p create_field just reads key given by
108  * @p field_name and creates instance using @p FieldBase<...>::function_factory.
109  * Function should return empty SharedField (that is shared_ptr to FieldBase).
110  *
111  * Implementation of these descendants is necessary:
112  * 1) for backward compatibility with old BCD input files
113  * 2) for setting pressure values are piezometric head values
114  */
115  /**
116  * Note for future:
117  * We pass through parameter @p field information about field that holds the factory which are necessary
118  * for interpreting user input and create particular field instance. It would be clearer to pass these information
119  * when the factory is assigned to a field. Moreover some information may not be set to field at all but directly passed
120  * to the factory.
121  */
122  class FactoryBase {
123  public:
124  /**
125  * Default method that creates an instance of FieldBase for field.
126  *
127  * Reads key given by @p field_name and creates the field instance using
128  * @p FieldBase<...>::function_factory.
129  */
130  virtual FieldBasePtr create_field(Input::Record rec, const FieldCommon &field);
131 
132  /**
133  * Check if Input::Record accessor contains data of field given by input_name.
134  *
135  * Returns true when ever the method create_field returns non-null pointer, otherwise returns false.
136  */
137  virtual bool is_active_field_descriptor(const Input::Record &in_rec, const std::string &input_name);
138  };
139 
140  /**
141  * Default constructor.
142  *
143  */
144  Field();
145 
146  Field(const string &name);
147 
148  /**
149  * Constructor that must be used for create of MultiField components.
150  *
151  * Set parameters @p component_index_, @p shared_->input_name_ and @p name_.
152  * Parameter name_ of Field is consisted of component name and MultiField name.
153  */
154  Field(unsigned int component_index, string input_name, string name = "");
155 
156  /**
157  * Copy constructor. Keeps shared history, declaration data, mesh.
158  */
159  Field(const Field &other);
160 
161  /**
162  * Assignment operator. Same properties as copy constructor, but class member name_ is not copied.
163  *
164  * Question: do we really need this, isn't copy constructor enough?
165  * Answer: It is necessary in (usual) case when Field instance is created as the class member
166  * but is filled later by assignment possibly from other class.
167  * TODO: operator can be merged with copy constructor, but we must provide to set correct value
168  * of name in method copy_from
169  */
170  Field &operator=(const Field &other);
171 
172 
173  /// Return appropriate value to BulkPoint in FieldValueCache
174  typename Value::return_type operator() (BulkPoint &p);
175 
176 
177  /// Return appropriate value to EdgePoint in FieldValueCache
178  typename Value::return_type operator() (SidePoint &p);
179 
180 
181  /**
182  * Returns reference to input type of particular field instance, this is static member @p input_type of the corresponding FieldBase class
183  * (i.e. with same template parameters). However, for fields returning "Enum" we have to create whole unique Input::Type hierarchy using following method
184  * @p make_input_tree.
185  * every instance since every such field use different Selection for initialization, even if all returns just unsigned int.
186  */
187  IT::Instance get_input_type() override;
188 
189  IT::Array get_multifield_input_type() override;
190 
191 
192  /**
193  * By this method you can allow that the field need not to be set on regions (and times) where the given @p control_field is
194  * FieldConstant and has value in given @p value_list. We check this in the set_time method. Through this mechanism we
195  * can switch of e.g. boundary data fields according to the type of the boundary condition.
196  */
197  auto disable_where(
198  const Field<spacedim, typename FieldValue<spacedim>::Enum > &control_field,
199  const vector<FieldEnum> &value_list) -> Field &;
200 
201 
202 
203  /**
204  * Set mesh pointer and resize region arrays.
205  *
206  * Implements abstract method.
207  */
208  void set_mesh(const Mesh &mesh) override;
209 
210 
211  /**
212  * Direct read access to the table of Field pointers on regions.
213  */
214  //std::shared_ptr< FieldBaseType > operator[] (Region reg);
215 
216  /**
217  * Implementation of @p FieldCommonBase::is_constant().
218  * See also Field<>::field_result which provide better information about special field values.
219  */
220  bool is_constant(Region reg) override;
221 
222  /**
223  * Assigns given @p field to all regions in region set given by @p region_set_names.
224  * Field is added to the history with given time and possibly used in the next call of the set_time method.
225  * Caller is responsible for correct construction of given field.
226  *
227  * Use this method only if necessary.
228  */
229  void set(FieldBasePtr field, double time, std::vector<std::string> region_set_names = {"ALL"});
230 
231  /**
232  * Same as before but the field is first created using FieldBase::function_factory(), from
233  * given abstract record accessor @p a_rec.
234  */
235  void set(const Input::AbstractRecord &a_rec, double time, std::vector<std::string> region_set_names = {"ALL"});
236 
237  /**
238  * Check that whole field list is set, possibly use default values for unset regions
239  * and call set_time for every field in the field list.
240  *
241  * Returns true if the field has been changed.
242  */
243  bool set_time(const TimeStep &time, LimitSide limit_side) override;
244 
245  /**
246  * Check that other has same type and assign from it.
247  */
248  void copy_from(const FieldCommon & other) override;
249 
250  /**
251  * Implementation of FieldCommonBase::output().
252  */
253  void field_output(std::shared_ptr<OutputTime> stream, OutputTime::DiscreteSpace type) override;
254 
255  /**
256  * Returns true, if field is currently set to a time in which it is discontinuous.
257  */
258  //bool is_jump_time();
259 
260  /**
261  * Check that the field is in fact FieldFE set on all bulk regions, return shared pointer to that FieldFE or NULL
262  * if the Field is not FieldFE.
263  */
264  std::shared_ptr< FieldFE<spacedim, Value> > get_field_fe();
265 
266 
267  /**
268  * @brief Indicates special field states.
269  *
270  * Return possible values from the enum @p FieldResult, see description there.
271  * The initial state is @p field_none, if the field is correctly set on all regions of the @p region_set given as parameter
272  * we return state @p field_other
273  * - Special field values spatially constant. Could allow optimization of tensor multiplication and
274  * tensor or vector addition. field_result_ should be set in constructor and in set_time method of particular Field implementation.
275  * We return value @p result_none, if the field is not initialized on the region of the given element accessor @p elm.
276  * Other possible results are: result_zeros, result_eye, result_ones, result_constant, result_other
277  * see @p FieldResult for explanation.
278  */
279  FieldResult field_result( RegionSet region_set) const override;
280 
281  /**
282  * Return value of input type attribute 'field_value_shape' that is appended to the
283  * input type of this field in FieldSet::make_field_descriptor_type and also to the output field selection
284  * created in EquationOutput::make_output_type.
285  * This attribute is used by GeoMop to have semantics of the input and output field data.
286  *
287  * Attribute value is a valid JSON (and/or flow style YAML) with keys:
288  * 'subfields' - value True for multifields, False or not presented for single value fields
289  * 'shape' - [ NRows, Ncols] ... given by FieldValue
290  * 'type' - <element type> (Double or Integer) ... given by FieldValue
291  * 'limit' - bounds of the field values.
292  *
293  */
294  std::string get_value_attribute() const override;
295 
296  /**
297  * Add a new factory for creating Field algorithms on individual regions.
298  * The last factory is tried first, the last one is always the default implementation
299  * Field<...>::FactoryBase.
300  *
301  * The Field<...> object keeps a list of such factories. When the instance of a new field algorithm
302  * has to be created from the input field descriptor, we pass through the list of factories backward
303  * and let factories to create the field algorithm instance from the actual input field descriptor.
304  * The first instance (non-null pointer) is used.
305  */
306  void add_factory(std::shared_ptr<FactoryBase> factory);
307 
308  void set_input_list(const Input::Array &list, const TimeGovernor &tg) override;
309 
310  /**
311  * Create and return shared_ptr to ElementDataCache appropriate to Field. Data cache is given by discrete @p space_type
312  * and is stored into data structures of output time @p stream for postponed output too.
313  */
314  void set_output_data_cache(OutputTime::DiscreteSpace space_type, std::shared_ptr<OutputTime> stream) override;
315 
316  /**
317  * Interpolate given field into output discrete @p space_type and store the values
318  * into storage of output time @p stream for postponed output.
319  */
320  void compute_field_data(OutputTime::DiscreteSpace space_type, std::shared_ptr<OutputTime> stream);
321 
322  /// Implements FieldCommon::cache_allocate
323  void cache_reallocate(const ElementCacheMap &cache_map, unsigned int region_idx) const override;
324 
325  /// Implements FieldCommon::cache_update
326  void cache_update(ElementCacheMap &cache_map, unsigned int region_patch_idx) const override;
327 
328  /// Implements FieldCommon::value_cache
329  FieldValueCache<double> * value_cache() override;
330 
331  /// Implements FieldCommon::value_cache
332  const FieldValueCache<double> * value_cache() const override;
333 
334  /**
335  * Implementation of FieldCommon::set_dependency().
336  */
337  std::vector<const FieldCommon *> set_dependency(unsigned int i_reg) const override;
338 
339  /// Implements FieldCommon::fill_data_value
340  void fill_data_value(const std::vector<int> &offsets) override;
341 
342  /// Implements FieldCommon::fill_observe_value
343  void fill_observe_value(std::shared_ptr<ElementDataCacheBase> output_cache_base, const std::vector<int> &offsets) override;
344 
345 protected:
346 
347  /// Return item of @p value_cache_ given by i_cache_point.
348  typename Value::return_type operator[] (unsigned int i_cache_point) const;
349 
350  /**
351  * Read input into @p regions_history_ possibly pop some old values from the
352  * history queue to keep its size less then @p history_length_limit_.
353  */
354  void update_history(const TimeStep &time);
355 
356  /**
357  * Check that whole field list (@p region_fields_) is set, possibly use default values for unset regions.
358  */
359  void check_initialized_region_fields_();
360 
361  /**************** Shared data **************/
362 
363  /// Pair: time, pointer to FieldBase instance
364  typedef pair<double, FieldBasePtr> HistoryPoint;
365  /// Nearest history of one region.
366  typedef boost::circular_buffer<HistoryPoint> RegionHistory;
367 
368  struct SharedData {
369 
370  /**
371  * History for every region. Shared among copies.
372  */
374  };
375 
376  /**************** Data per copy **************/
377 
378  std::shared_ptr<SharedData> data_;
379 
380  /**
381  * If this pointer is set, turn off check of initialization in the
382  * @p set_time method on the regions where the method @p get_constant_enum_value
383  * of the control field returns value from @p no_check_values_. This
384  * field is private copy, its set_time method is called from the
385  * set_Time method of actual object.
386  */
388  std::shared_ptr<ControlField> no_check_control_field_;
389 
390  /**
391  * Table with pointers to fields on individual regions.
392  */
394 
396 
397  /**
398  * Field value data cache
399  *
400  * Data is ordered like three dimensional table. The highest level is determinated by subsets,
401  * those data ranges are holds in subset_starts. Data block size of each subset is determined
402  * by number of eval_points (of subset) and maximal number of stored elements.
403  * The table is allocated to hold all subsets, but only those marked in used_subsets are updated.
404  * Order of subsets is same as in eval_points.
405  */
407 
408  /// ElementDataCache used during field output, object is shared with OutputTime
409  std::shared_ptr<ElementDataCache<typename Value::element_type>> output_data_cache_;
410 
411 
412 
413  template<int dim, class Val>
414  friend class MultiField;
415 
416  template< typename CALLABLE, typename TUPLE, int INDEX >
418 
419 };
420 
421 
422 
423 
424 
425 
426 
427 #endif /* FIELD_HH_ */
Definitions of ASSERTS.
Base point accessor class.
Definition: eval_subset.hh:55
Directing class of FieldValueCache.
Class holds local coordinations of evaluating points (bulk and sides) specified by element dimension.
Definition: eval_points.hh:43
Space< spacedim >::Point Point
Common abstract parent of all Field<...> classes.
Definition: field_common.hh:77
Container for various descendants of FieldCommonBase.
Definition: field_set.hh:159
Class template representing a field with values dependent on: point, element, and region.
Definition: field.hh:92
FieldValueCache< typename Value::element_type > value_cache_
Definition: field.hh:406
FieldAlgorithmBase< spacedim, Value >::Point Point
Definition: field.hh:97
Field< spacedim, typename FieldValue< spacedim >::Enum > ControlField
Definition: field.hh:387
FieldAlgorithmBase< spacedim, Value > FieldBaseType
Definition: field.hh:95
pair< double, FieldBasePtr > HistoryPoint
Pair: time, pointer to FieldBase instance.
Definition: field.hh:364
boost::circular_buffer< HistoryPoint > RegionHistory
Nearest history of one region.
Definition: field.hh:366
Value ValueType
Definition: field.hh:98
std::shared_ptr< ElementDataCache< typename Value::element_type > > output_data_cache_
ElementDataCache used during field output, object is shared with OutputTime.
Definition: field.hh:409
std::shared_ptr< SharedData > data_
Definition: field.hh:378
std::shared_ptr< ControlField > no_check_control_field_
Definition: field.hh:388
std::vector< std::shared_ptr< FactoryBase > > factories_
Definition: field.hh:395
std::shared_ptr< FieldBaseType > FieldBasePtr
Definition: field.hh:96
std::vector< FieldBasePtr > region_fields_
Definition: field.hh:393
Accessor to the polymorphic input data of a type given by an AbstracRecord object.
Definition: accessors.hh:458
Accessor to input data conforming to declared Array.
Definition: accessors.hh:566
Accessor to the data with type Type::Record.
Definition: accessors.hh:291
Class for declaration of inputs sequences.
Definition: type_base.hh:339
Helper class that stores data of generic types.
Definition: type_generic.hh:89
Definition: mesh.h:362
Class for representation of a vector of fields of the same physical quantity.
Definition: multi_field.hh:87
General point a+ side_begin_ + ccessor allow iterate over quadrature points of given side defined in ...
Definition: eval_subset.hh:116
Basic time management functionality for unsteady (and steady) solvers (class Equation).
Representation of one time step..
FieldResult
LimitSide
Definition: field_common.hh:64
@ Value
Definition: field.hh:62
std::vector< RegionHistory > region_history_
Definition: field.hh:373
base case for building up arguments for the function call
Definition: field_model.hh:64
Basic time management class.