GetFEM  5.4.2
getfem_models.h
Go to the documentation of this file.
1 /* -*- c++ -*- (enables emacs c++ mode) */
2 /*===========================================================================
3 
4  Copyright (C) 2009-2020 Yves Renard
5 
6  This file is a part of GetFEM
7 
8  GetFEM is free software; you can redistribute it and/or modify it
9  under the terms of the GNU Lesser General Public License as published
10  by the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version along with the GCC Runtime Library
12  Exception either version 3.1 or (at your option) any later version.
13  This program is distributed in the hope that it will be useful, but
14  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  License and GCC Runtime Library Exception for more details.
17  You should have received a copy of the GNU Lesser General Public License
18  along with this program; if not, write to the Free Software Foundation,
19  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
20 
21  As a special exception, you may use this file as it is a part of a free
22  software library without restriction. Specifically, if other files
23  instantiate templates or use macros or inline functions from this file,
24  or you compile this file and link it with other files to produce an
25  executable, this file does not by itself cause the resulting executable
26  to be covered by the GNU Lesser General Public License. This exception
27  does not however invalidate any other reasons why the executable file
28  might be covered by the GNU Lesser General Public License.
29 
30 ===========================================================================*/
31 
32 /**
33  @file getfem_models.h
34  @author Yves Renard <[email protected]>
35  @date March 21, 2009.
36  @brief Model representation in Getfem.
37 */
38 
39 #ifndef GETFEM_MODELS_H__
40 #define GETFEM_MODELS_H__
41 
42 #include "getfem_generic_assembly.h"
43 #include "getfem_assembling.h"
44 #include "getfem_partial_mesh_fem.h"
45 #include "getfem_im_data.h"
46 
47 namespace getfem {
48 
49  class virtual_brick;
50  /** type of pointer on a brick */
51  typedef std::shared_ptr<const virtual_brick> pbrick;
52 
53  class virtual_dispatcher;
54  typedef std::shared_ptr<const virtual_dispatcher> pdispatcher;
55 
56  class virtual_time_scheme;
57  typedef std::shared_ptr<const virtual_time_scheme> ptime_scheme;
58 
59  // Event management : The model has to react when something has changed in
60  // the context and ask for corresponding (linear) bricks to recompute
61  // some terms.
62  // For the moment two events are taken into account
63  // - Change in a mesh_fem
64  // - Change in the data of a variable
65  // For this, a brick has to declare on which variable it depends and
66  // on which data. When a linear brick depend on a variable, the
67  // recomputation is done when the eventual corresponding mesh_fem
68  // is changed (or the size of the variable for a fixed size variable).
69  // When a linear brick depend on a data, the recomputation is done
70  // when the corresponding vector value is changed. If a variable is used
71  // as a data, it has to be declared as a data by the brick.
72  // A nonlinear brick is recomputed at each assembly of the tangent system.
73  // Remember this behavior when some changed are done on the variable
74  // and/or data.
75  // The change on a mesh_im is not taken into account for the moment.
76  // The different versions of the variables is not taken into account
77  // separately.
78 
79 
80 
81 
82  //=========================================================================
83  //
84  // Model object.
85  //
86  //=========================================================================
87 
88  // For backward compatibility with version 3.0
89  typedef model_real_plain_vector modeling_standard_plain_vector;
90  typedef model_real_sparse_vector modeling_standard_sparse_vector;
91  typedef model_real_sparse_matrix modeling_standard_sparse_matrix;
92  typedef model_complex_plain_vector modeling_standard_complex_plain_vector;
93  typedef model_complex_sparse_vector modeling_standard_complex_sparse_vector;
94  typedef model_complex_sparse_matrix modeling_standard_complex_sparse_matrix;
95 
96 
97  /** A prefix to refer to the previous version of a variable*/
98  const auto PREFIX_OLD = std::string{"Old_"};
99  const auto PREFIX_OLD_LENGTH = 4;
100 
101  /** Does the variable have Old_ prefix*/
102  bool is_old(const std::string &name);
103 
104  /** Strip the variable name from prefix Old_ if it has one*/
105  std::string no_old_prefix_name(const std::string &name);
106 
107  std::string sup_previous_and_dot_to_varname(std::string v);
108 
109  /** ``Model'' variables store the variables, the data and the
110  description of a model. This includes the global tangent matrix, the
111  right hand side and the constraints. There are two kinds of models, the
112  ``real'' and the ``complex'' models.
113  */
114  class APIDECL model : public context_dependencies,
115  virtual public dal::static_stored_object {
116 
117  protected:
118 
119  // State variables of the model
120  bool complex_version;
121  bool is_linear_;
122  bool is_symmetric_;
123  bool is_coercive_;
124  mutable model_real_sparse_matrix
125  rTM, // tangent matrix (only primary variables), real version
126  internal_rTM; // coupling matrix between internal and primary vars (no empty rows)
127  mutable model_complex_sparse_matrix cTM; // tangent matrix, complex version
128  mutable model_real_plain_vector
129  rrhs, // residual vector of primary variables (after condensation)
130  full_rrhs, // residual vector of primary and internal variables (pre-condensation)
131  internal_sol; // partial solution for internal variables (after condensation)
132  mutable model_complex_plain_vector crhs;
133  mutable bool act_size_to_be_done;
134  dim_type leading_dim;
135  getfem::lock_factory locks_;
136 
137  // Variables and parameters of the model
138 
139  enum var_description_filter {
140  VDESCRFILTER_NO = 0, // Variable being directly the dofs of a given fem
141  VDESCRFILTER_REGION = 1, /* Variable being the dofs of a fem on a mesh
142  * region (uses mf.dof_on_region). */
143  VDESCRFILTER_INFSUP = 2, /* Variable being the dofs of a fem on a mesh
144  * region with an additional filter on a mass
145  * matrix with respect to another fem. */
146  VDESCRFILTER_CTERM = 4, /* Variable being the dofs of a fem on a mesh
147  * region with an additional filter with the
148  * coupling term with respect to another variable.*/
149  VDESCRFILTER_REGION_CTERM = 5, /* both region and cterm. */
150  }; // INFSUP and CTERM are incompatible
151 
152  struct var_description {
153 
154  bool is_variable; // This is a variable or a parameter.
155  bool is_disabled; // For a variable, to be solved or not
156  bool is_complex; // The variable is complex numbers
157  bool is_affine_dependent; // The variable depends in an affine way
158  // to another variable.
159  bool is_internal; // An internal variable defined on integration
160  // points, condensed out of the global system.
161  bool is_fem_dofs; // The variable is the dofs of a fem
162  size_type n_iter; // Number of versions of the variable stored.
163  size_type n_temp_iter; // Number of additional temporary versions
164  size_type default_iter; // default iteration number.
165 
166  ptime_scheme ptsc; // For optional time integration scheme
167 
168  var_description_filter filter; // Version of an optional filter
169  // on the dofs
170  size_type filter_region; // Optional region for the filter
171  std::string filter_var; // Optional variable name for the filter
172  mesh_im const *filter_mim; // Optional integration method for the filter
173 
174  // fem or im_data description of the variable
175  mesh_fem const *mf; // Main fem of the variable.
176  ppartial_mesh_fem partial_mf; // Filter with respect to mf.
177  im_data const *imd; // im data description
178 
179  bgeot::multi_index qdims; // For data having a qdim != of the fem
180  // (dim per dof for dof data)
181  // and for constant variables.
182  gmm::uint64_type v_num;
183  std::vector<gmm::uint64_type> v_num_data;
184 
185  gmm::sub_interval I; // For a variable : indices on the whole system.
186  // For an affine dependent variable, should be the same as the
187  // original variable
188  std::vector<model_real_plain_vector> real_value;
189  std::vector<model_complex_plain_vector> complex_value;
190  std::vector<gmm::uint64_type> v_num_var_iter;
191  std::vector<gmm::uint64_type> v_num_iter;
192 
193  // For affine dependent variables
194  model_real_plain_vector affine_real_value;
195  model_complex_plain_vector affine_complex_value;
196  scalar_type alpha; // Factor for the affine dependent variables
197  std::string org_name; // Name of the original variable for affine
198  // dependent variables
199 
200  var_description(bool is_var = false, bool is_compl = false,
201  mesh_fem const *mf_ = 0, im_data const *imd_ = 0,
202  size_type n_it = 1,
203  var_description_filter filter_ = VDESCRFILTER_NO,
204  size_type filter_reg = size_type(-1),
205  const std::string &filter_var_ = std::string(""),
206  mesh_im const *filter_mim_ = 0)
207  : is_variable(is_var), is_disabled(false), is_complex(is_compl),
208  is_affine_dependent(false), is_internal(false),
209  is_fem_dofs(mf_ != 0),
210  n_iter(std::max(size_type(1), n_it)), n_temp_iter(0),
211  default_iter(0), ptsc(0),
212  filter(filter_), filter_region(filter_reg), filter_var(filter_var_),
213  filter_mim(filter_mim_), mf(mf_), imd(imd_), qdims(),
214  v_num(0), v_num_data(n_iter, act_counter()), I(0,0), alpha(1)
215  {
216  if (filter != VDESCRFILTER_NO && mf != 0)
217  partial_mf = std::make_shared<partial_mesh_fem>(*mf);
218  // v_num_data = v_num;
219  if (qdims.size() == 0) qdims.push_back(1);
220  GMM_ASSERT1(qdim(), "Attempt to create a null size variable");
221  }
222 
223  size_type qdim() const { return qdims.total_size(); }
224 
225  // add a temporary version for time integration schemes. Automatically
226  // set the default iter to it. id_num is an identifier. Do not add
227  // the version if a temporary already exist with this identifier.
228  size_type add_temporary(gmm::uint64_type id_num);
229 
230  void clear_temporaries();
231 
232  const mesh_fem &associated_mf() const {
233  GMM_ASSERT1(is_fem_dofs, "This variable is not linked to a fem");
234  return (filter == VDESCRFILTER_NO) ? *mf : *partial_mf;
235  }
236 
237  const mesh_fem *passociated_mf() const {
238  if (is_fem_dofs)
239  return (filter == VDESCRFILTER_NO || partial_mf.get() == 0)
240  ? mf : partial_mf.get();
241  return 0;
242  }
243 
244  size_type size() const { // Should control that the variable is
245  // indeed initialized by actualize_sizes() ...
246  return is_complex ? complex_value[0].size()
247  : real_value[0].size();
248  }
249  inline bool is_enabled() const { return !is_disabled; }
250 
251  void set_size();
252  }; // struct var_description
253 
254  public:
255 
256  typedef std::vector<std::string> varnamelist;
257  typedef std::vector<const mesh_im *> mimlist;
258  typedef std::vector<model_real_sparse_matrix> real_matlist;
259  typedef std::vector<model_complex_sparse_matrix> complex_matlist;
260  typedef std::vector<model_real_plain_vector> real_veclist;
261  typedef std::vector<model_complex_plain_vector> complex_veclist;
262 
263  struct term_description {
264  bool is_matrix_term; // tangent matrix term or rhs term.
265  bool is_symmetric; // Term have to be symmetrized.
266  bool is_global; // Specific global term for highly coupling bricks
267  std::string var1, var2;
268 
269  term_description(const std::string &v)
270  : is_matrix_term(false), is_symmetric(false),
271  is_global(false), var1(sup_previous_and_dot_to_varname(v)) {}
272  term_description(const std::string &v1, const std::string &v2,
273  bool issym)
274  : is_matrix_term(true), is_symmetric(issym), is_global(false),
275  var1(sup_previous_and_dot_to_varname(v1)), var2(v2) {}
276  term_description(bool ism, bool issym)
277  : is_matrix_term(ism), is_symmetric(issym), is_global(true) {}
278  };
279 
280  typedef std::vector<term_description> termlist;
281 
282  enum build_version {
283  BUILD_RHS = 1,
284  BUILD_MATRIX = 2,
285  BUILD_ALL = 3,
286  BUILD_ON_DATA_CHANGE = 4,
287  BUILD_WITH_LIN = 8, // forced calculation of linear terms
288  BUILD_RHS_WITH_LIN = 9, // = BUILD_RHS | BUILD_WITH_LIN_RHS
289  BUILD_WITH_INTERNAL = 16,
290  BUILD_RHS_WITH_INTERNAL = 17, // = BUILD_RHS | BUILD_WITH_INTERNAL
291  BUILD_MATRIX_CONDENSED = 18, // = BUILD_MATRIX | BUILD_WITH_INTERNAL
292  BUILD_ALL_CONDENSED = 19, // = BUILD_ALL | BUILD_WITH_INTERNAL
293  };
294 
295  protected:
296 
297  // rmatlist and cmatlist could had been csc_matrix vectors to reduce the
298  // amount of memory (but this would require an additional copy).
299  struct brick_description {
300  mutable bool terms_to_be_computed;
301  mutable gmm::uint64_type v_num;
302  pbrick pbr; // brick pointer
303  pdispatcher pdispatch; // Optional dispatcher
304  size_type nbrhs; // Additional rhs for dispatcher.
305  varnamelist vlist; // List of variables used by the brick.
306  varnamelist dlist; // List of data used by the brick.
307  termlist tlist; // List of terms build by the brick
308  mimlist mims; // List of integration methods.
309  size_type region; // Optional region size_type(-1) for all.
310  bool is_update_brick; // Flag for declaring special type of brick
311  // with no contributions.
312  mutable scalar_type external_load; // External load computed in assembly
313 
314  mutable model_real_plain_vector coeffs;
315  mutable scalar_type matrix_coeff = scalar_type(0);
316  // Matrices the brick has to fill in (real version)
317  mutable real_matlist rmatlist;
318  // Rhs the brick has to fill in (real version)
319  mutable std::vector<real_veclist> rveclist;
320  // additional rhs for symmetric terms (real version)
321  mutable std::vector<real_veclist> rveclist_sym;
322  // Matrices the brick has to fill in (complex version)
323  mutable complex_matlist cmatlist;
324  // Rhs the brick has to fill in (complex version)
325  mutable std::vector<complex_veclist> cveclist;
326  // additional rhs for symmetric terms (complex version)
327  mutable std::vector<complex_veclist> cveclist_sym;
328 
329  brick_description() : v_num(0) {}
330 
331  brick_description(pbrick p, const varnamelist &vl,
332  const varnamelist &dl, const termlist &tl,
333  const mimlist &mms, size_type reg)
334  : terms_to_be_computed(true), v_num(0), pbr(p), pdispatch(0), nbrhs(1),
335  vlist(vl), dlist(dl), tlist(tl), mims(mms), region(reg),
336  is_update_brick(false), external_load(0),
337  rveclist(1), rveclist_sym(1), cveclist(1),
338  cveclist_sym(1) { }
339  };
340 
341  typedef std::map<std::string, var_description> VAR_SET;
342  mutable VAR_SET variables; // Variables list of the model
343  std::vector<brick_description> bricks; // Bricks list of the model
344  dal::bit_vector valid_bricks, active_bricks;
345  std::map<std::string, pinterpolate_transformation> transformations;
346  std::map<std::string, pelementary_transformation> elem_transformations;
347  std::map<std::string, psecondary_domain> secondary_domains;
348 
349  // Structure dealing with time integration scheme
350  int time_integration; // 0 : no, 1 : time step, 2 : init
351  bool init_step;
352  scalar_type time_step; // Time step (dt) for time integration schemes
353  scalar_type init_time_step; // Time step for initialization of derivatives
354 
355  // Structure dealing with simple dof constraints
356  typedef std::map<size_type, scalar_type> real_dof_constraints_var;
357  typedef std::map<size_type, complex_type> complex_dof_constraints_var;
358  mutable std::map<std::string, real_dof_constraints_var>
359  real_dof_constraints;
360  mutable std::map<std::string, complex_dof_constraints_var>
361  complex_dof_constraints;
362  void clear_dof_constraints()
363  { real_dof_constraints.clear(); complex_dof_constraints.clear(); }
364 
365  // Structure dealing with nonlinear expressions
366  struct gen_expr {
367  std::string expr;
368  const mesh_im &mim;
369  size_type region;
370  std::string secondary_domain;
371  gen_expr(const std::string &expr_, const mesh_im &mim_,
372  size_type region_, const std::string &secdom)
373  : expr(expr_), mim(mim_), region(region_), secondary_domain(secdom) {}
374  };
375 
376  // Structure for assignment in assembly
377  struct assignement_desc {
378  std::string varname;
379  std::string expr;
380  size_type region;
381  bool before;
382  size_type order;
383  };
384 
385  std::list<assignement_desc> assignments;
386 
387  mutable std::list<gen_expr> generic_expressions;
388 
389  // Groups of variables for interpolation on different meshes
390  // generic assembly
391  std::map<std::string, std::vector<std::string> > variable_groups;
392 
393  ga_macro_dictionary macro_dict;
394 
395 
396  virtual void actualize_sizes() const;
397  bool check_name_validity(const std::string &name, bool assert=true) const;
398  void brick_init(size_type ib, build_version version,
399  size_type rhs_ind = 0) const;
400 
401  void init() { complex_version = false; act_size_to_be_done = false; }
402 
403  void resize_global_system() const;
404 
405  //to be performed after to_variables is called.
406  virtual void post_to_variables_step();
407 
408  scalar_type approx_external_load_; // Computed by assembly procedure
409  // with BUILD_RHS option.
410 
411  VAR_SET::const_iterator find_variable(const std::string &name) const;
412  const var_description &variable_description(const std::string &name) const;
413 
414  public:
415 
416  void add_generic_expression(const std::string &expr, const mesh_im &mim,
417  size_type region,
418  const std::string &secondary_domain = "") const {
419  generic_expressions.push_back(gen_expr(expr, mim, region,
420  secondary_domain));
421  }
422  void add_external_load(size_type ib, scalar_type e) const
423  { bricks[ib].external_load = e; }
424  scalar_type approx_external_load() { return approx_external_load_; }
425  // call the brick if necessary
426  void update_brick(size_type ib, build_version version) const;
427  void linear_brick_add_to_rhs(size_type ib, size_type ind_data,
428  size_type n_iter) const;
429  void update_affine_dependent_variables();
430  void brick_call(size_type ib, build_version version,
431  size_type rhs_ind = 0) const;
432  model_real_plain_vector &rhs_coeffs_of_brick(size_type ib) const
433  { return bricks[ib].coeffs; }
434  scalar_type &matrix_coeff_of_brick(size_type ib) const
435  { return bricks[ib].matrix_coeff; }
436  bool is_var_newer_than_brick(const std::string &varname,
437  size_type ib, size_type niter = size_type(-1)) const;
438  bool is_var_mf_newer_than_brick(const std::string &varname,
439  size_type ib) const;
440  bool is_mim_newer_than_brick(const mesh_im &mim,
441  size_type ib) const;
442 
443  pbrick brick_pointer(size_type ib) const {
444  GMM_ASSERT1(valid_bricks[ib], "Inexistent brick");
445  return bricks[ib].pbr;
446  }
447 
448  void variable_list(varnamelist &vl) const
449  { for (const auto &v : variables) vl.push_back(v.first); }
450 
451  void define_variable_group(const std::string &group_name,
452  const std::vector<std::string> &nl);
453  bool variable_group_exists(const std::string &group_name) const
454  { return variable_groups.count(group_name) > 0; }
455 
456  const std::vector<std::string> &
457  variable_group(const std::string &group_name) const {
458  GMM_ASSERT1(variable_group_exists(group_name),
459  "Undefined variable group " << group_name);
460  return (variable_groups.find(group_name))->second;
461  }
462 
463  void clear_assembly_assignments(void) { assignments.clear(); }
464  void add_assembly_assignments(const std::string &dataname,
465  const std::string &expr,
466  size_type rg = size_type(-1),
467  size_type order = 1,
468  bool before = false);
469 
470  /* function to be called by Dirichlet bricks */
471  void add_real_dof_constraint(const std::string &varname, size_type dof,
472  scalar_type val) const
473  { real_dof_constraints[varname][dof] = val; }
474  /* function to be called by Dirichlet bricks */
475  void add_complex_dof_constraint(const std::string &varname, size_type dof,
476  complex_type val) const
477  { complex_dof_constraints[varname][dof] = val; }
478 
479 
480  void add_temporaries(const varnamelist &vl, gmm::uint64_type id_num) const;
481 
482  const mimlist &mimlist_of_brick(size_type ib) const
483  { return bricks[ib].mims; }
484 
485  const varnamelist &varnamelist_of_brick(size_type ib) const
486  { return bricks[ib].vlist; }
487 
488  const varnamelist &datanamelist_of_brick(size_type ib) const
489  { return bricks[ib].dlist; }
490 
491  size_type region_of_brick(size_type ib) const
492  { return bricks[ib].region; }
493 
494  bool temporary_uptodate(const std::string &varname,
495  gmm::uint64_type id_num, size_type &ind) const;
496 
497  size_type n_iter_of_variable(const std::string &name) const {
498  return variables.count(name) == 0 ? size_type(0)
499  : variables[name].n_iter;
500  }
501 
502  void set_default_iter_of_variable(const std::string &varname,
503  size_type ind) const;
504  void reset_default_iter_of_variables(const varnamelist &vl) const;
505 
506  void update_from_context() const { act_size_to_be_done = true; }
507 
508  const model_real_sparse_matrix &linear_real_matrix_term
509  (size_type ib, size_type iterm);
510 
511  const model_complex_sparse_matrix &linear_complex_matrix_term
512  (size_type ib, size_type iterm);
513 
514  /** Disable a brick. */
515  void disable_brick(size_type ib) {
516  GMM_ASSERT1(valid_bricks[ib], "Inexistent brick");
517  active_bricks.del(ib);
518  }
519 
520  /** Enable a brick. */
521  void enable_brick(size_type ib) {
522  GMM_ASSERT1(valid_bricks[ib], "Inexistent brick");
523  active_bricks.add(ib);
524  }
525 
526  /** Disable a variable (and its attached mutlipliers). */
527  void disable_variable(const std::string &name);
528 
529  /** Enable a variable (and its attached mutlipliers). */
530  void enable_variable(const std::string &name, bool enabled=true);
531 
532  /** States if a name corresponds to a declared variable. */
533  bool variable_exists(const std::string &name) const;
534 
535  /** States if a variable is disabled (treated as data). */
536  bool is_disabled_variable(const std::string &name) const;
537 
538  /** States if a name corresponds to a declared data or disabled variable. */
539  bool is_data(const std::string &name) const;
540 
541  /** States if a name corresponds to a declared data. */
542  bool is_true_data(const std::string &name) const;
543 
544  /** States if a variable is condensed out of the global system. */
545  bool is_internal_variable(const std::string &name) const;
546 
547  bool is_affine_dependent_variable(const std::string &name) const;
548 
549  const std::string &org_variable(const std::string &name) const;
550 
551  const scalar_type &factor_of_variable(const std::string &name) const;
552 
553  void set_factor_of_variable(const std::string &name, scalar_type a);
554 
555  bool is_im_data(const std::string &name) const;
556 
557  const im_data *pim_data_of_variable(const std::string &name) const;
558 
559  const gmm::uint64_type &
560  version_number_of_data_variable(const std::string &varname,
561  size_type niter = size_type(-1)) const;
562 
563  /** Boolean which says if the model deals with real or complex unknowns
564  and data. */
565  bool is_complex() const { return complex_version; }
566 
567  /** Return true if all the model terms do not affect the coercivity of
568  the whole tangent system. */
569  bool is_coercive() const { return is_coercive_; }
570 
571  /** Return true if the model has at least one internal variable. */
572  bool has_internal_variables() const {
573  for (const auto &v : variables)
574  if (v.second.is_internal && !v.second.is_disabled) return true;
575  return false;
576  }
577 
578  /** Return true if all the model terms do not affect the coercivity of
579  the whole tangent system. */
580  bool is_symmetric() const { return is_symmetric_; }
581 
582  /** Return true if all the model terms are linear. */
583  bool is_linear() const { return is_linear_; }
584 
585  /** Total number of degrees of freedom in the model. */
586  size_type nb_dof(bool with_internal=false) const;
587 
588  /** Number of internal degrees of freedom in the model. */
589  size_type nb_internal_dof() const { return nb_dof(true)-nb_dof(); }
590 
591  /** Number of primary degrees of freedom in the model. */
592  size_type nb_primary_dof() const { return nb_dof(); }
593 
594  /** Leading dimension of the meshes used in the model. */
595  dim_type leading_dimension() const { return leading_dim; }
596 
597  /** Gives a non already existing variable name begining by `name`. */
598  std::string new_name(const std::string &name);
599 
600  const gmm::sub_interval &
601  interval_of_variable(const std::string &name) const;
602 
603  /** Gives the access to the vector value of a variable. For the real
604  version. */
605  const model_real_plain_vector &
606  real_variable(const std::string &name, size_type niter) const;
607 
608  /**The same as above, but either accessing the latest variable version,
609  or the previous, if using "Old_" prefix*/
610  const model_real_plain_vector &
611  real_variable(const std::string &name) const;
612 
613  /** Gives the access to the vector value of a variable. For the complex
614  version. */
615  const model_complex_plain_vector &
616  complex_variable(const std::string &name, size_type niter) const;
617 
618  /**The same as above, but either accessing the latest variable version,
619  or the previous, if using "Old_" prefix*/
620  const model_complex_plain_vector &
621  complex_variable(const std::string &name) const;
622 
623  /** Gives the write access to the vector value of a variable. Make a
624  change flag of the variable set. For the real version. */
625  model_real_plain_vector &
626  set_real_variable(const std::string &name, size_type niter) const;
627 
628  /**The same as above, but for either latest variable, or
629  for the previous, if prefixed with "Old_".*/
630  model_real_plain_vector &
631  set_real_variable(const std::string &name) const;
632 
633  /** Gives the write access to the vector value of a variable. Make a
634  change flag of the variable set. For the complex version. */
635  model_complex_plain_vector &
636  set_complex_variable(const std::string &name, size_type niter) const;
637 
638  /**The same as above, but either accessing the latest variable version,
639  or the previous, if using "Old_" prefix*/
640  model_complex_plain_vector &
641  set_complex_variable(const std::string &name) const;
642 
643  model_real_plain_vector &
644  set_real_constant_part(const std::string &name) const;
645 
646  model_complex_plain_vector &
647  set_complex_constant_part(const std::string &name) const;
648 
649  private:
650  template<typename VECTOR, typename T>
651  void from_variables(VECTOR &V, bool with_internal, T) const {
652  for (const auto &v : variables)
653  if (v.second.is_variable && !v.second.is_affine_dependent
654  && !v.second.is_disabled
655  && (with_internal || !v.second.is_internal))
656  gmm::copy(v.second.real_value[0], gmm::sub_vector(V, v.second.I));
657  }
658 
659  template<typename VECTOR, typename T>
660  void from_variables(VECTOR &V, bool with_internal, std::complex<T>) const {
661  for (const auto &v : variables)
662  if (v.second.is_variable && !v.second.is_affine_dependent
663  && !v.second.is_disabled
664  && (with_internal || !v.second.is_internal))
665  gmm::copy(v.second.complex_value[0], gmm::sub_vector(V, v.second.I));
666  }
667 
668  public:
669  template<typename VECTOR>
670  void from_variables(VECTOR &V, bool with_internal=false) const {
671  typedef typename gmm::linalg_traits<VECTOR>::value_type T;
672  context_check(); if (act_size_to_be_done) actualize_sizes();
673  from_variables(V, with_internal, T());
674  }
675 
676  private:
677  template<typename VECTOR, typename T>
678  void to_variables(const VECTOR &V, bool with_internal, T) {
679  for (auto &&v : variables)
680  if (v.second.is_variable && !v.second.is_affine_dependent
681  && !v.second.is_disabled
682  && (with_internal || !v.second.is_internal)) {
683  gmm::copy(gmm::sub_vector(V, v.second.I), v.second.real_value[0]);
684  v.second.v_num_data[0] = act_counter();
685  }
686  update_affine_dependent_variables();
687  this->post_to_variables_step();
688  }
689 
690  template<typename VECTOR, typename T>
691  void to_variables(const VECTOR &V, bool with_internal, std::complex<T>) {
692  for (auto &&v : variables)
693  if (v.second.is_variable && !v.second.is_affine_dependent
694  && !v.second.is_disabled
695  && (with_internal || !v.second.is_internal)) {
696  gmm::copy(gmm::sub_vector(V, v.second.I), v.second.complex_value[0]);
697  v.second.v_num_data[0] = act_counter();
698  }
699  update_affine_dependent_variables();
700  this->post_to_variables_step();
701  }
702 
703  public:
704  template<typename VECTOR>
705  void to_variables(const VECTOR &V, bool with_internal=false) {
706  typedef typename gmm::linalg_traits<VECTOR>::value_type T;
707  context_check(); if (act_size_to_be_done) actualize_sizes();
708  to_variables(V, with_internal, T());
709  }
710 
711  /** Add a fixed size variable to the model assumed to be a vector.
712  niter is the number of version of the variable stored. */
713  void add_fixed_size_variable(const std::string &name, size_type size,
714  size_type niter = 1);
715 
716  /** Add a fixed size variable to the model whith given tensor dimensions.
717  niter is the number of version of the variable stored. */
718  void add_fixed_size_variable(const std::string &name,
719  const bgeot::multi_index &sizes,
720  size_type niter = 1);
721 
722  /** Add a fixed size data to the model. niter is the number of version
723  of the data stored, for time integration schemes. */
724  void add_fixed_size_data(const std::string &name, size_type size,
725  size_type niter = 1);
726 
727  /** Add a fixed size data to the model. niter is the number of version
728  of the data stored, for time integration schemes. */
729  void add_fixed_size_data(const std::string &name,
730  const bgeot::multi_index &sizes,
731  size_type niter = 1);
732 
733  /** Resize a fixed size variable (or data) of the model. */
734  void resize_fixed_size_variable(const std::string &name, size_type size);
735 
736  /** Resize a fixed size variable (or data) of the model. */
737  void resize_fixed_size_variable(const std::string &name,
738  const bgeot::multi_index &sizes);
739 
740  /** Add a fixed size data (assumed to be a vector) to the model and
741  initialized with v. */
742  template <typename VECT>
743  void add_initialized_fixed_size_data(const std::string &name,
744  const VECT &v) {
745  this->add_fixed_size_data(name, gmm::vect_size(v));
746  if (this->is_complex())
747  gmm::copy(v, this->set_complex_variable(name));
748  else
749  gmm::copy(gmm::real_part(v), this->set_real_variable(name));
750  }
751 
752  /** Add a fixed size data (assumed to be a vector) to the model and
753  initialized with v. */
754  template <typename VECT>
755  void add_initialized_fixed_size_data(const std::string &name,
756  const VECT &v,
757  const bgeot::multi_index &sizes) {
758  this->add_fixed_size_data(name, sizes);
759  if (this->is_complex())
760  gmm::copy(v, this->set_complex_variable(name));
761  else
762  gmm::copy(gmm::real_part(v), this->set_real_variable(name));
763  }
764 
765  /** Add a fixed size data (assumed to be a matrix) to the model and
766  initialized with M. */
767  void add_initialized_matrix_data(const std::string &name,
768  const base_matrix &M);
769  void add_initialized_matrix_data(const std::string &name,
770  const base_complex_matrix &M);
771 
772  /** Add a fixed size data (assumed to be a tensor) to the model and
773  initialized with t. */
774  void add_initialized_tensor_data(const std::string &name,
775  const base_tensor &t);
776  void add_initialized_tensor_data(const std::string &name,
777  const base_complex_tensor &t);
778 
779 
780  /** Add a scalar data (i.e. of size 1) to the model initialized with e. */
781  template <typename T>
782  void add_initialized_scalar_data(const std::string &name, T e) {
783  this->add_fixed_size_data(name, 1, 1);
784  if (this->is_complex())
785  this->set_complex_variable(name)[0] = e;
786  else
787  this->set_real_variable(name)[0] = gmm::real(e);
788  }
789 
790 
791  /** Add variable defined at integration points. */
792  void add_im_variable(const std::string &name, const im_data &imd,
793  size_type niter = 1);
794  /** Add internal variable, defined at integration points and condensed. */
795  void add_internal_im_variable(const std::string &name,
796  const im_data &imd);
797  /** Add data defined at integration points. */
798  void add_im_data(const std::string &name, const im_data &imd,
799  size_type niter = 1);
800 
801  /** Add a variable being the dofs of a finite element method to the model.
802  niter is the number of version of the variable stored, for time
803  integration schemes. */
804  void add_fem_variable(const std::string &name, const mesh_fem &mf,
805  size_type niter = 1);
806 
807  /** Add a variable linked to a fem with the dof filtered with respect
808  to a mesh region. Only the dof returned by the dof_on_region
809  method of `mf` will be kept. niter is the number of version
810  of the data stored, for time integration schemes. */
811  void add_filtered_fem_variable(const std::string &name, const mesh_fem &mf,
812  size_type region, size_type niter = 1);
813 
814 
815  /** Add a "virtual" variable be an affine depedent variable with respect
816  to another variable. Mainly used for time integration scheme for
817  instance to represent time derivative of variables.
818  `alpha` is the multiplicative scalar of the dependency. */
819  void add_affine_dependent_variable(const std::string &name,
820  const std::string &org_name,
821  scalar_type alpha = scalar_type(1));
822 
823  /** Add a data being the dofs of a finite element method to the model.*/
824  void add_fem_data(const std::string &name, const mesh_fem &mf,
825  dim_type qdim = 1, size_type niter = 1);
826 
827  /** Add a data being the dofs of a finite element method to the model.*/
828  void add_fem_data(const std::string &name, const mesh_fem &mf,
829  const bgeot::multi_index &sizes, size_type niter = 1);
830 
831  /** Add an initialized fixed size data to the model, assumed to be a
832  vector field if the size of the vector is a multiple of the dof
833  number. */
834  template <typename VECT>
835  void add_initialized_fem_data(const std::string &name, const mesh_fem &mf,
836  const VECT &v) {
837  this->add_fem_data(name, mf,
838  dim_type(gmm::vect_size(v) / mf.nb_dof()), 1);
839  if (this->is_complex())
840  gmm::copy(v, this->set_complex_variable(name));
841  else
842  gmm::copy(gmm::real_part(v), this->set_real_variable(name));
843  }
844 
845  /** Add a fixed size data to the model. The data is a tensor of given
846  sizes on each dof of the finite element method. */
847  template <typename VECT>
848  void add_initialized_fem_data(const std::string &name, const mesh_fem &mf,
849  const VECT &v,
850  const bgeot::multi_index &sizes) {
851  this->add_fem_data(name, mf, sizes, 1);
852  if (this->is_complex())
853  gmm::copy(v, this->set_complex_variable(name));
854  else
855  gmm::copy(gmm::real_part(v), this->set_real_variable(name));
856  }
857 
858  /** Add a particular variable linked to a fem being a multiplier with
859  respect to a primal variable. The dof will be filtered with the
860  gmm::range_basis function applied on the terms of the model which
861  link the multiplier and the primal variable. Optimized for boundary
862  multipliers. niter is the number of version of the data stored,
863  for time integration schemes. */
864  void add_multiplier(const std::string &name, const mesh_fem &mf,
865  const std::string &primal_name,
866  size_type niter = 1);
867 
868  /** Add a particular variable linked to a fem being a multiplier with
869  respect to a primal variable and a region. The dof will be filtered
870  both with the gmm::range_basis function applied on the terms of
871  the model which link the multiplier and the primal variable and on
872  the dof on the given region. Optimized for boundary
873  multipliers. niter is the number of version of the data stored,
874  for time integration schemes. */
875  void add_multiplier(const std::string &name, const mesh_fem &mf,
876  size_type region, const std::string &primal_name,
877  size_type niter = 1);
878 
879  /** Add a particular variable linked to a fem being a multiplier with
880  respect to a primal variable. The dof will be filtered with the
881  gmm::range_basis function applied on the mass matrix between the fem
882  of the multiplier and the one of the primal variable.
883  Optimized for boundary multipliers. niter is the number of version
884  of the data stored, for time integration schemes. */
885  void add_multiplier(const std::string &name, const mesh_fem &mf,
886  const std::string &primal_name, const mesh_im &mim,
887  size_type region, size_type niter = 1);
888 
889  /** Dictonnary of user defined macros. */
890  const ga_macro_dictionary &macro_dictionary() const { return macro_dict; }
891 
892  /** Add a macro definition for the high generic assembly language.
893  This macro can be used for the definition of generic assembly bricks.
894  The name of a macro cannot coincide with a variable name. */
895  void add_macro(const std::string &name, const std::string &expr);
896 
897  /** Delete a previously defined macro definition. */
898  void del_macro(const std::string &name);
899 
900  /** Says if a macro of that name has been defined. */
901  bool macro_exists(const std::string &name) const
902  { return macro_dict.macro_exists(name); }
903 
904  /** Delete a variable or data of the model. */
905  void delete_variable(const std::string &varname);
906 
907  /** Gives the access to the mesh_fem of a variable if any. Throw an
908  exception otherwise. */
909  const mesh_fem &mesh_fem_of_variable(const std::string &name) const;
910 
911  /** Gives a pointer to the mesh_fem of a variable if any. 0 otherwise.*/
912  const mesh_fem *pmesh_fem_of_variable(const std::string &name) const;
913 
914 
915  bgeot::multi_index qdims_of_variable(const std::string &name) const;
916  size_type qdim_of_variable(const std::string &name) const;
917 
918  /** Gives the access to the tangent matrix. For the real version. */
919  const model_real_sparse_matrix &
920  real_tangent_matrix(bool internal=false) const {
921  GMM_ASSERT1(!complex_version, "This model is a complex one");
922  context_check(); if (act_size_to_be_done) actualize_sizes();
923  return internal ? internal_rTM : rTM;
924  }
925 
926  /** Gives the access to the tangent matrix. For the complex version. */
927  const model_complex_sparse_matrix &complex_tangent_matrix() const {
928  GMM_ASSERT1(complex_version, "This model is a real one");
929  context_check(); if (act_size_to_be_done) actualize_sizes();
930  return cTM;
931  }
932 
933  /** Gives access to the right hand side of the tangent linear system.
934  For the real version. An assembly of the rhs has to be done first. */
935  const model_real_plain_vector &real_rhs(bool with_internal=false) const {
936  GMM_ASSERT1(!complex_version, "This model is a complex one");
937  context_check(); if (act_size_to_be_done) actualize_sizes();
938  return (with_internal && gmm::vect_size(full_rrhs)) ? full_rrhs : rrhs;
939  }
940 
941  /** Gives write access to the right hand side of the tangent linear system.
942  Some solvers need to manipulate the model rhs directly so that for
943  example internal condensed variables can be treated properly. */
944  model_real_plain_vector &set_real_rhs(bool with_internal=false) const {
945  GMM_ASSERT1(!complex_version, "This model is a complex one");
946  context_check(); if (act_size_to_be_done) actualize_sizes();
947  return (with_internal && gmm::vect_size(full_rrhs)) ? full_rrhs : rrhs;
948  }
949 
950  /** Gives access to the partial solution for condensed internal variables.
951  A matrix assembly with condensation has to be done first. */
952  const model_real_plain_vector &internal_solution() const {
953  GMM_ASSERT1(!complex_version, "This model is a complex one");
954  context_check(); if (act_size_to_be_done) actualize_sizes();
955  return internal_sol;
956  }
957 
958  /** Gives access to the part of the right hand side of a term of
959  a particular nonlinear brick. Does not account of the eventual time
960  dispatcher. An assembly of the rhs has to be done first.
961  For the real version. */
962  const model_real_plain_vector &real_brick_term_rhs
963  (size_type ib, size_type ind_term = 0, bool sym = false,
964  size_type ind_iter = 0) const
965  {
966  GMM_ASSERT1(!complex_version, "This model is a complex one");
967  context_check(); if (act_size_to_be_done) actualize_sizes();
968  GMM_ASSERT1(valid_bricks[ib], "Inexistent brick");
969  GMM_ASSERT1(ind_term < bricks[ib].tlist.size(), "Inexistent term");
970  GMM_ASSERT1(ind_iter < bricks[ib].nbrhs, "Inexistent iter");
971  GMM_ASSERT1(!sym || bricks[ib].tlist[ind_term].is_symmetric,
972  "Term is not symmetric");
973  if (sym)
974  return bricks[ib].rveclist_sym[ind_iter][ind_term];
975  else
976  return bricks[ib].rveclist[ind_iter][ind_term];
977  }
978 
979  /** Gives access to the right hand side of the tangent linear system.
980  For the complex version. */
981  const model_complex_plain_vector &complex_rhs() const {
982  GMM_ASSERT1(complex_version, "This model is a real one");
983  context_check(); if (act_size_to_be_done) actualize_sizes();
984  return crhs;
985  }
986 
987  /** Gives write access to the right hand side of the tangent linear system.
988  Some solvers need to manipulate the model rhs directly so that for
989  example internal condensed variables can be treated properly. */
990  model_complex_plain_vector &set_complex_rhs() const {
991  GMM_ASSERT1(complex_version, "This model is a real one");
992  context_check(); if (act_size_to_be_done) actualize_sizes();
993  return crhs;
994  }
995 
996  /** Gives access to the part of the right hand side of a term of a
997  particular nonlinear brick. Does not account of the eventual time
998  dispatcher. An assembly of the rhs has to be done first.
999  For the complex version. */
1000  const model_complex_plain_vector &complex_brick_term_rhs
1001  (size_type ib, size_type ind_term = 0, bool sym = false,
1002  size_type ind_iter = 0) const
1003  {
1004  GMM_ASSERT1(!complex_version, "This model is a complex one");
1005  context_check(); if (act_size_to_be_done) actualize_sizes();
1006  GMM_ASSERT1(valid_bricks[ib], "Inexistent brick");
1007  GMM_ASSERT1(ind_term < bricks[ib].tlist.size(), "Inexistent term");
1008  GMM_ASSERT1(ind_iter < bricks[ib].nbrhs, "Inexistent iter");
1009  GMM_ASSERT1(!sym || bricks[ib].tlist[ind_term].is_symmetric,
1010  "Term is not symmetric");
1011  if (sym)
1012  return bricks[ib].cveclist_sym[ind_iter][ind_term];
1013  else
1014  return bricks[ib].cveclist[ind_iter][ind_term];
1015  }
1016 
1017  /** List the model variables and constant. */
1018  void listvar(std::ostream &ost) const;
1019 
1020  void listresiduals(std::ostream &ost) const;
1021 
1022  /** List the model bricks. */
1023  void listbricks(std::ostream &ost, size_type base_id = 0) const;
1024 
1025  /** Return the model brick ids. */
1026  const dal::bit_vector& get_active_bricks() const {
1027  return active_bricks;
1028  }
1029 
1030  /** Force the re-computation of a brick for the next assembly. */
1031  void touch_brick(size_type ib) {
1032  GMM_ASSERT1(valid_bricks[ib], "Inexistent brick");
1033  bricks[ib].terms_to_be_computed = true;
1034  }
1035 
1036  /** Add a brick to the model. varname is the list of variable used
1037  and datanames the data used. If a variable is used as a data, it
1038  should be declared in the datanames (it will depend on the value of
1039  the variable not only on the fem). Returns the brick index. */
1040  size_type add_brick(pbrick pbr, const varnamelist &varnames,
1041  const varnamelist &datanames,
1042  const termlist &terms, const mimlist &mims,
1043  size_type region);
1044 
1045  /** Delete the brick of index ib from the model. */
1046  void delete_brick(size_type ib);
1047 
1048  /** Add an integration method to a brick. */
1049  void add_mim_to_brick(size_type ib, const mesh_im &mim);
1050 
1051  /** Change the term list of a brick. Used for very special bricks only. */
1052  void change_terms_of_brick(size_type ib, const termlist &terms);
1053 
1054  /** Change the variable list of a brick. Used for very special bricks only.
1055  */
1056  void change_variables_of_brick(size_type ib, const varnamelist &vl);
1057 
1058  /** Change the data list of a brick. Used for very special bricks only.
1059  */
1060  void change_data_of_brick(size_type ib, const varnamelist &vl);
1061 
1062  /** Change the mim list of a brick. Used for very special bricks only.
1063  */
1064  void change_mims_of_brick(size_type ib, const mimlist &ml);
1065 
1066  /** Change the update flag of a brick. Used for very special bricks only.
1067  */
1068  void change_update_flag_of_brick(size_type ib, bool flag);
1069 
1070  void set_time(scalar_type t = scalar_type(0), bool to_init = true);
1071 
1072  scalar_type get_time();
1073 
1074  void set_time_step(scalar_type dt) { time_step = dt; }
1075  scalar_type get_time_step() const { return time_step; }
1076  scalar_type get_init_time_step() const { return init_time_step; }
1077  int is_time_integration() const { return time_integration; }
1078  void set_time_integration(int ti) { time_integration = ti; }
1079  bool is_init_step() const { return init_step; }
1080  void cancel_init_step() { init_step = false; }
1081  void call_init_affine_dependent_variables(int version);
1082  void shift_variables_for_time_integration();
1083  void copy_init_time_derivative();
1084  void add_time_integration_scheme(const std::string &varname,
1085  ptime_scheme ptsc);
1086  void perform_init_time_derivative(scalar_type ddt)
1087  { init_step = true; init_time_step = ddt; }
1088 
1089 
1090  /** Add a time dispacther to a brick. */
1091  void add_time_dispatcher(size_type ibrick, pdispatcher pdispatch);
1092 
1093  void set_dispatch_coeff();
1094 
1095  /** For transient problems. Initialisation of iterations. */
1096  virtual void first_iter();
1097 
1098  /** For transient problems. Prepare the next iterations. In particular
1099  shift the version of the variables.
1100  */
1101  virtual void next_iter();
1102 
1103  /** Add an interpolate transformation to the model to be used with the
1104  generic assembly.
1105  */
1106  void add_interpolate_transformation(const std::string &name,
1107  pinterpolate_transformation ptrans) {
1108  if (secondary_domain_exists(name))
1109  GMM_ASSERT1(false, "An secondary domain with the same "
1110  "name already exists");
1111  if (transformations.count(name) > 0)
1112  GMM_ASSERT1(name.compare("neighbor_element"), "neighbor_element is a "
1113  "reserved interpolate transformation name");
1114  transformations[name] = ptrans;
1115  }
1116 
1117  /** Get a pointer to the interpolate transformation `name`.
1118  */
1119  pinterpolate_transformation
1120  interpolate_transformation(const std::string &name) const {
1121  std::map<std::string, pinterpolate_transformation>::const_iterator
1122  it = transformations.find(name);
1123  GMM_ASSERT1(it != transformations.end(), "Inexistent transformation " << name);
1124  return it->second;
1125  }
1126 
1127  /** Tests if `name` corresponds to an interpolate transformation.
1128  */
1129  bool interpolate_transformation_exists(const std::string &name) const
1130  { return transformations.count(name) > 0; }
1131 
1132  /** Add an elementary transformation to the model to be used with the
1133  generic assembly.
1134  */
1135  void add_elementary_transformation(const std::string &name,
1136  pelementary_transformation ptrans) {
1137  elem_transformations[name] = ptrans;
1138  }
1139 
1140  /** Get a pointer to the elementary transformation `name`.
1141  */
1142  pelementary_transformation
1143  elementary_transformation(const std::string &name) const {
1144  std::map<std::string, pelementary_transformation>::const_iterator
1145  it = elem_transformations.find(name);
1146  GMM_ASSERT1(it != elem_transformations.end(),
1147  "Inexistent elementary transformation " << name);
1148  return it->second;
1149  }
1150 
1151  /** Tests if `name` corresponds to an elementary transformation.
1152  */
1153  bool elementary_transformation_exists(const std::string &name) const
1154  { return elem_transformations.count(name) > 0; }
1155 
1156 
1157  /** Add a secondary domain to the model to be used with the
1158  generic assembly.
1159  */
1160  void add_secondary_domain(const std::string &name,
1161  psecondary_domain ptrans) {
1162  if (interpolate_transformation_exists(name))
1163  GMM_ASSERT1(false, "An interpolate transformation with the same "
1164  "name already exists");secondary_domains[name] = ptrans;
1165  }
1166 
1167  /** Get a pointer to the interpolate transformation `name`.
1168  */
1169  psecondary_domain
1170  secondary_domain(const std::string &name) const {
1171  auto it = secondary_domains.find(name);
1172  GMM_ASSERT1(it != secondary_domains.end(),
1173  "Inexistent transformation " << name);
1174  return it->second;
1175  }
1176 
1177  /** Tests if `name` corresponds to an interpolate transformation.
1178  */
1179  bool secondary_domain_exists(const std::string &name) const
1180  { return secondary_domains.count(name) > 0; }
1181 
1182  /** Gives the name of the variable of index `ind_var` of the brick
1183  of index `ind_brick`. */
1184  const std::string &varname_of_brick(size_type ind_brick,
1185  size_type ind_var);
1186 
1187  /** Gives the name of the data of index `ind_data` of the brick
1188  of index `ind_brick`. */
1189  const std::string &dataname_of_brick(size_type ind_brick,
1190  size_type ind_data);
1191 
1192  /** Assembly of the tangent system taking into account the terms
1193  from all bricks. */
1194  virtual void assembly(build_version version);
1195 
1196  /** Gives the assembly string corresponding to the Neumann term of
1197  the fem variable `varname` on `region`. It is deduced from the
1198  assembly string declared by the model bricks.
1199  `region` should be the index of a boundary region
1200  on the mesh where `varname` is defined. Care to call this function
1201  only after all the volumic bricks have been declared.
1202  Complains, if a brick omit to declare an assembly string.
1203  */
1204  std::string Neumann_term(const std::string &varname, size_type region);
1205 
1206  virtual void clear();
1207 
1208  explicit model(bool comp_version = false);
1209 
1210  /** check consistency of RHS and Stiffness matrix for brick with
1211  @param ind_brick - index of the brick
1212  */
1213  void check_brick_stiffness_rhs(size_type ind_brick) const;
1214 
1215 
1216  };
1217 
1218  //=========================================================================
1219  //
1220  // Time integration scheme object.
1221  //
1222  //=========================================================================
1223 
1224  /** The time integration scheme object provides the necessary methods
1225  for the model object to apply a time integration scheme to an
1226  evolutionnary problem.
1227  **/
1228  class APIDECL virtual_time_scheme {
1229 
1230  public:
1231 
1232  virtual void init_affine_dependent_variables(model &md) const = 0;
1233  virtual void init_affine_dependent_variables_precomputation(model &md)
1234  const = 0;
1235  virtual void time_derivative_to_be_initialized
1236  (std::string &name_v, std::string &name_previous_v) const = 0;
1237  virtual void shift_variables(model &md) const = 0;
1238  virtual ~virtual_time_scheme() {}
1239  };
1240 
1241  void add_theta_method_for_first_order(model &md, const std::string &varname,
1242  scalar_type theta);
1243 
1244  void add_theta_method_for_second_order(model &md, const std::string &varname,
1245  scalar_type theta);
1246 
1247  void add_Newmark_scheme(model &md, const std::string &varname,
1248  scalar_type beta, scalar_type gamma);
1249 
1250  void add_Houbolt_scheme(model &md, const std::string &varname);
1251 
1252 
1253 
1254  //=========================================================================
1255  //
1256  // Time dispatcher object.
1257  //
1258  //=========================================================================
1259 
1260  /** The time dispatcher object modify the result of a brick in order to
1261  apply a time integration scheme.
1262  **/
1263  class APIDECL virtual_dispatcher {
1264 
1265  protected:
1266 
1267  size_type nbrhs_;
1268  std::vector<std::string> param_names;
1269 
1270  public:
1271 
1272  size_type nbrhs() const { return nbrhs_; }
1273 
1274  typedef model::build_version build_version;
1275 
1276  virtual void set_dispatch_coeff(const model &, size_type) const
1277  { GMM_ASSERT1(false, "Time dispatcher with not set_dispatch_coeff !"); }
1278 
1279  virtual void next_real_iter
1280  (const model &, size_type, const model::varnamelist &,
1281  const model::varnamelist &,
1282  model::real_matlist &,
1283  std::vector<model::real_veclist> &,
1284  std::vector<model::real_veclist> &,
1285  bool) const {
1286  GMM_ASSERT1(false, "Time dispatcher with not defined first real iter !");
1287  }
1288 
1289  virtual void next_complex_iter
1290  (const model &, size_type, const model::varnamelist &,
1291  const model::varnamelist &,
1292  model::complex_matlist &,
1293  std::vector<model::complex_veclist> &,
1294  std::vector<model::complex_veclist> &,
1295  bool) const{
1296  GMM_ASSERT1(false,"Time dispatcher with not defined first comples iter");
1297  }
1298 
1299  virtual void asm_real_tangent_terms
1300  (const model &, size_type,
1301  model::real_matlist &, std::vector<model::real_veclist> &,
1302  std::vector<model::real_veclist> &,
1303  build_version) const {
1304  GMM_ASSERT1(false, "Time dispatcher with not defined real tangent "
1305  "terms !");
1306  }
1307 
1308  virtual void asm_complex_tangent_terms
1309  (const model &, size_type,
1310  model::complex_matlist &, std::vector<model::complex_veclist> &,
1311  std::vector<model::complex_veclist> &,
1312  build_version) const {
1313  GMM_ASSERT1(false, "Time dispatcher with not defined complex tangent "
1314  "terms !");
1315  }
1316 
1317  virtual_dispatcher(size_type _nbrhs) : nbrhs_(_nbrhs) {
1318  GMM_ASSERT1(_nbrhs > 0, "Time dispatcher with no rhs");
1319  }
1320  virtual ~virtual_dispatcher() {}
1321 
1322  };
1323 
1324  // ----------------------------------------------------------------------
1325  //
1326  // theta-method dispatcher
1327  //
1328  // ----------------------------------------------------------------------
1329 
1330  class APIDECL theta_method_dispatcher : public virtual_dispatcher {
1331 
1332  public:
1333 
1334  typedef model::build_version build_version;
1335 
1336  void set_dispatch_coeff(const model &md, size_type ib) const;
1337 
1338  template <typename MATLIST, typename VECTLIST>
1339  inline void next_iter(const model &md, size_type ib,
1340  const model::varnamelist &/* vl */,
1341  const model::varnamelist &/* dl */,
1342  MATLIST &/* matl */,
1343  VECTLIST &vectl, VECTLIST &vectl_sym,
1344  bool first_iter) const {
1345  if (first_iter) md.update_brick(ib, model::BUILD_RHS);
1346 
1347  // shift the rhs
1348  for (size_type i = 0; i < vectl[0].size(); ++i)
1349  gmm::copy(vectl[0][i], vectl[1][i]);
1350  for (size_type i = 0; i < vectl_sym[0].size(); ++i)
1351  gmm::copy(vectl_sym[0][i], vectl_sym[1][i]);
1352 
1353  // add the component represented by the linear matrix terms to the
1354  // supplementary rhs
1355  md.linear_brick_add_to_rhs(ib, 1, 0);
1356  }
1357 
1358  void next_real_iter
1359  (const model &md, size_type ib, const model::varnamelist &vl,
1360  const model::varnamelist &dl, model::real_matlist &matl,
1361  std::vector<model::real_veclist> &vectl,
1362  std::vector<model::real_veclist> &vectl_sym, bool first_iter) const;
1363 
1364  void next_complex_iter
1365  (const model &md, size_type ib, const model::varnamelist &vl,
1366  const model::varnamelist &dl,
1367  model::complex_matlist &matl,
1368  std::vector<model::complex_veclist> &vectl,
1369  std::vector<model::complex_veclist> &vectl_sym,
1370  bool first_iter) const;
1371 
1372  void asm_real_tangent_terms
1373  (const model &md, size_type ib, model::real_matlist &/* matl */,
1374  std::vector<model::real_veclist> &/* vectl */,
1375  std::vector<model::real_veclist> &/* vectl_sym */,
1376  build_version version) const;
1377 
1378  virtual void asm_complex_tangent_terms
1379  (const model &md, size_type ib, model::complex_matlist &/* matl */,
1380  std::vector<model::complex_veclist> &/* vectl */,
1381  std::vector<model::complex_veclist> &/* vectl_sym */,
1382  build_version version) const;
1383 
1384  theta_method_dispatcher(const std::string &THETA);
1385  };
1386 
1387  //=========================================================================
1388  //
1389  // Functions adding standard time dispatchers.
1390  //
1391  //=========================================================================
1392 
1393  /** Add a theta-method time dispatcher to a list of bricks. For instance,
1394  a matrix term $K$ will be replaced by
1395  $\theta K U^{n+1} + (1-\theta) K U^{n}$.
1396  */
1397  void APIDECL add_theta_method_dispatcher(model &md, dal::bit_vector ibricks,
1398  const std::string &THETA);
1399 
1400  /** Function which udpate the velocity $v^{n+1}$ after the computation
1401  of the displacement $u^{n+1}$ and before the next iteration. Specific
1402  for theta-method and when the velocity is included in the data
1403  of the model.
1404  */
1405  void APIDECL velocity_update_for_order_two_theta_method
1406  (model &md, const std::string &U, const std::string &V,
1407  const std::string &pdt, const std::string &ptheta);
1408 
1409 
1410  /** Add a midpoint time dispatcher to a list of bricks. For instance,
1411  a nonlinear term $K(U)$ will be replaced by
1412  $K((U^{n+1} + U^{n})/2)$.
1413  */
1414  void APIDECL add_midpoint_dispatcher(model &md, dal::bit_vector ibricks);
1415 
1416  /** Function which udpate the velocity $v^{n+1}$ after the computation
1417  of the displacement $u^{n+1}$ and before the next iteration. Specific
1418  for Newmark scheme and when the velocity is included in the data
1419  of the model. This version inverts the mass matrix by a conjugate
1420  gradient.
1421  */
1422  void APIDECL velocity_update_for_Newmark_scheme
1423  (model &md, size_type id2dt2b, const std::string &U, const std::string &V,
1424  const std::string &pdt, const std::string &ptwobeta,
1425  const std::string &pgamma);
1426 
1427  //=========================================================================
1428  //
1429  // Brick object.
1430  //
1431  //=========================================================================
1432 
1433  /** The virtual brick has to be derived to describe real model bricks.
1434  The set_flags method has to be called by the derived class.
1435  The virtual methods asm_real_tangent_terms and/or
1436  asm_complex_tangent_terms have to be defined.
1437  The brick should not store data. The data have to be stored in the
1438  model object.
1439  **/
1440  class APIDECL virtual_brick {
1441  protected:
1442  bool islinear; // The brick add a linear term or not.
1443  bool issymmetric; // The brick add a symmetric term or not.
1444  bool iscoercive; // The brick add a potentially coercive terms or not.
1445  // (in particular, not a term involving a multiplier)
1446  bool isreal; // The brick admits a real version or not.
1447  bool iscomplex; // The brick admits a complex version or not.
1448  bool isinit; // internal flag.
1449  bool compute_each_time; // The brick is linear but needs to be computed
1450  // each time it is evaluated.
1451  bool isUpdateBrick; // The brick does not contribute any terms to the
1452  // system matrix or right-hand side, but only updates state variables.
1453  std::string name; // Name of the brick.
1454 
1455  public:
1456 
1457  typedef model::build_version build_version;
1458 
1459  virtual_brick() { isinit = false; }
1460  virtual ~virtual_brick() { }
1461  void set_flags(const std::string &bname, bool islin, bool issym,
1462  bool iscoer, bool ire, bool isco, bool each_time = false) {
1463  name = bname;
1464  islinear = islin; issymmetric = issym; iscoercive = iscoer;
1465  isreal = ire; iscomplex = isco; isinit = true;
1466  compute_each_time = each_time;
1467  }
1468 
1469 # define BRICK_NOT_INIT GMM_ASSERT1(isinit, "Set brick flags !")
1470  bool is_linear() const { BRICK_NOT_INIT; return islinear; }
1471  bool is_symmetric() const { BRICK_NOT_INIT; return issymmetric; }
1472  bool is_coercive() const { BRICK_NOT_INIT; return iscoercive; }
1473  bool is_real() const { BRICK_NOT_INIT; return isreal; }
1474  bool is_complex() const { BRICK_NOT_INIT; return iscomplex; }
1475  bool is_to_be_computed_each_time() const
1476  { BRICK_NOT_INIT; return compute_each_time; }
1477  const std::string &brick_name() const { BRICK_NOT_INIT; return name; }
1478 
1479 
1480  /** Assembly of bricks real tangent terms.
1481  In case of Getfem's compilation with OpenMP option,
1482  this method is executed on multiple threads.
1483  The parallelism is provided by distributing all tangent matrices and
1484  vectors and accumulating them later into the original. Additionally,
1485  by default, all mesh_region objects, participating in the assembly,
1486  are also partitioned. In order to avoid data race conditions, this
1487  method should not modify any data simultaneously accessible from
1488  multiple threads. In case this is unavoidable, the race can be
1489  prevented by distributing this data (of type T) between the threads
1490  via getfem::omp_distribute<T> (prefered method) or
1491  protected from concurrent access with mutexes (e.g. getfem::omp_lock)
1492  or OpenMP critical section. */
1493  virtual void asm_real_tangent_terms(const model &, size_type,
1494  const model::varnamelist &,
1495  const model::varnamelist &,
1496  const model::mimlist &,
1497  model::real_matlist &,
1498  model::real_veclist &,
1499  model::real_veclist &,
1500  size_type, build_version) const
1501  { /** doesn't have to be overriden if serial pre- post- assemblies are
1502  defined */
1503  }
1504 
1505 
1506  /** Assembly of bricks complex tangent terms.
1507  In case of Getfem's compilation with OpenMP option,
1508  this method is executed on multiple threads.
1509  The parallelism is provided by distributing all tangent matrices and
1510  vectors and accumulating them later into the original. Additionally,
1511  by default, all mesh_region objects, participating in the assembly,
1512  are also partitioned. In order to avoid data race conditions, this
1513  method should not modify any data simultaneously accessible from
1514  multiple threads. In case this is unavoidable, the race can be
1515  prevented by distributing this data (of type T) between the threads
1516  via getfem::omp_distribute<T> (prefered method) or
1517  protected from concurrent access with mutexes (e.g. getfem::omp_lock)
1518  or OpenMP critical section. */
1519  virtual void asm_complex_tangent_terms(const model &, size_type,
1520  const model::varnamelist &,
1521  const model::varnamelist &,
1522  const model::mimlist &,
1523  model::complex_matlist &,
1524  model::complex_veclist &,
1525  model::complex_veclist &,
1526  size_type, build_version) const
1527  { /** doesn't have to be overriden if serial pre- post- assemblies are
1528  defined*/
1529  }
1530 
1531 
1532  /** Peform any pre assembly action for real term assembly. The purpose of
1533  this method is to do any action that cannot be peformed in the main
1534  assembly routines in parallel.
1535  Possible action can be modification of the model object, cashing
1536  some data that cannot be distributed etc. */
1537  virtual void real_pre_assembly_in_serial(const model &, size_type,
1538  const model::varnamelist &,
1539  const model::varnamelist &,
1540  const model::mimlist &,
1541  model::real_matlist &,
1542  model::real_veclist &,
1543  model::real_veclist &,
1544  size_type, build_version) const { };
1545 
1546  /** Peform any pre assembly action for complex term assembly. The purpose
1547  of this method is to do any action that cannot be peformed in the
1548  main assembly routines in parallel.
1549  Possible action can be modification of the model object, cashing
1550  some data that cannot be distributed etc. */
1551  virtual void complex_pre_assembly_in_serial(const model &, size_type,
1552  const model::varnamelist &,
1553  const model::varnamelist &,
1554  const model::mimlist &,
1555  model::complex_matlist &,
1556  model::complex_veclist &,
1557  model::complex_veclist &,
1558  size_type, build_version) const { };
1559 
1560  /** Peform any post assembly action for real terms. The purpose of this
1561  method is to do any action that cannot be peformed in the main
1562  assembly routines in parallel.
1563  Possible action can be modification of the model object, cashing
1564  some data that cannot be distributed etc. */
1565  virtual void real_post_assembly_in_serial(const model &, size_type,
1566  const model::varnamelist &,
1567  const model::varnamelist &,
1568  const model::mimlist &,
1569  model::real_matlist &,
1570  model::real_veclist &,
1571  model::real_veclist &,
1572  size_type, build_version) const { };
1573 
1574  /** Peform any post assembly action for complex terms. The purpose of this
1575  method is to do any action that cannot be peformed in the main
1576  assembly routines in parallel.
1577  Possible action can be modification of the model object, cashing
1578  some data that cannot be distributed etc. */
1579  virtual void complex_post_assembly_in_serial(const model &, size_type,
1580  const model::varnamelist &,
1581  const model::varnamelist &,
1582  const model::mimlist &,
1583  model::complex_matlist &,
1584  model::complex_veclist &,
1585  model::complex_veclist &,
1586  size_type, build_version) const { };
1587 
1588 
1589  /** check consistency of stiffness matrix and rhs */
1590  void check_stiffness_matrix_and_rhs(const model &, size_type,
1591  const model::termlist& tlist,
1592  const model::varnamelist &,
1593  const model::varnamelist &,
1594  const model::mimlist &,
1595  model::real_matlist &,
1596  model::real_veclist &,
1597  model::real_veclist &, size_type rg,
1598  const scalar_type delta = 1e-8) const;
1599  /** The brick may declare an assembly string for the computation of the
1600  Neumann terms (in order to prescribe boundary conditions with
1601  Nitche's method). */
1602  virtual std::string declare_volume_assembly_string
1603  (const model &, size_type, const model::varnamelist &,
1604  const model::varnamelist &) const {
1605  GMM_ASSERT1(false, "No assemby string declared, computation of Neumann "
1606  "term impossible for brick " << name);
1607  }
1608 
1609  private:
1610  /** simultaneous call to real_pre_assembly, real_assembly
1611  and real_post_assembly */
1612  void full_asm_real_tangent_terms_(const model &, size_type,
1613  const model::varnamelist &,
1614  const model::varnamelist &,
1615  const model::mimlist &,
1616  model::real_matlist &,
1617  model::real_veclist &,
1618  model::real_veclist &,
1619  size_type, build_version) const;
1620  };
1621 
1622  //=========================================================================
1623  //
1624  // Functions adding standard bricks to the model.
1625  //
1626  //=========================================================================
1627 
1628  /** Add a term given by the weak form language expression `expr` which will
1629  be assembled in region `region` and with the integration method `mim`.
1630  Only the matrix term will be taken into account, assuming that it is
1631  linear.
1632  The advantage of declaring a term linear instead of nonlinear is that
1633  it will be assembled only once and no assembly is necessary for the
1634  residual.
1635  Take care that if the expression contains some variables and if the
1636  expression is a potential or of first order (i.e. describe the weak
1637  form, not the derivative of the weak form), the expression will be
1638  derivated with respect to all variables.
1639  You can specify if the term is symmetric, coercive or not.
1640  If you are not sure, the better is to declare the term not symmetric
1641  and not coercive. But some solvers (conjugate gradient for instance)
1642  are not allowed for non-coercive problems.
1643  `brickname` is an optional name for the brick.
1644  */
1645  size_type APIDECL add_linear_term
1646  (model &md, const mesh_im &mim, const std::string &expr,
1647  size_type region = size_type(-1), bool is_sym = false,
1648  bool is_coercive = false, const std::string &brickname = "",
1649  bool return_if_nonlin = false);
1650 
1651  inline size_type APIDECL add_linear_generic_assembly_brick
1652  (model &md, const mesh_im &mim, const std::string &expr,
1653  size_type region = size_type(-1), bool is_sym = false,
1654  bool is_coercive = false, const std::string &brickname = "",
1655  bool return_if_nonlin = false) {
1656  return add_linear_term(md, mim, expr, region, is_sym,
1657  is_coercive, brickname, return_if_nonlin);
1658  }
1659 
1660  /** Add a nonlinear term given by the weak form language expression `expr`
1661  which will be assembled in region `region` and with the integration
1662  method `mim`.
1663  The expression can describe a potential or a weak form. Second order
1664  terms (i.e. containing second order test functions, Test2) are not
1665  allowed.
1666  You can specify if the term is symmetric, coercive or not.
1667  If you are not sure, the better is to declare the term not symmetric
1668  and not coercive. But some solvers (conjugate gradient for instance)
1669  are not allowed for non-coercive problems.
1670  `brickname` is an optional name for the brick.
1671  */
1672  size_type APIDECL add_nonlinear_term
1673  (model &md, const mesh_im &mim, const std::string &expr,
1674  size_type region = size_type(-1), bool is_sym = false,
1675  bool is_coercive = false, const std::string &brickname = "");
1676 
1677  inline size_type APIDECL add_nonlinear_generic_assembly_brick
1678  (model &md, const mesh_im &mim, const std::string &expr,
1679  size_type region = size_type(-1), bool is_sym = false,
1680  bool is_coercive = false, const std::string &brickname = "") {
1681  return add_nonlinear_term(md, mim, expr, region,
1682  is_sym, is_coercive, brickname);
1683  }
1684 
1685 
1686  /** Add a source term given by the assembly string `expr` which will
1687  be assembled in region `region` and with the integration method `mim`.
1688  Only the residual term will be taken into account.
1689  Take care that if the expression contains some variables and if the
1690  expression is a potential, the expression will be
1691  derivated with respect to all variables.
1692  `brickname` is an optional name for the brick.
1693  */
1694  size_type APIDECL add_source_term
1695  (model &md, const mesh_im &mim, const std::string &expr,
1696  size_type region = size_type(-1),
1697  const std::string &brickname = std::string(),
1698  const std::string &directvarname = std::string(),
1699  const std::string &directdataname = std::string(),
1700  bool return_if_nonlin = false);
1701  inline size_type APIDECL add_source_term_generic_assembly_brick
1702  (model &md, const mesh_im &mim, const std::string &expr,
1703  size_type region = size_type(-1),
1704  const std::string &brickname = std::string(),
1705  const std::string &directvarname = std::string(),
1706  const std::string &directdataname = std::string(),
1707  bool return_if_nonlin = false) {
1708  return add_source_term(md, mim, expr, region, brickname,
1709  directvarname, directdataname, return_if_nonlin);
1710  }
1711 
1712  /** Adds a linear term given by a weak form language expression like
1713  ``add_linear_term`` function but for an integration on a direct
1714  product of two domains, a first specfied by ``mim`` and ``region``
1715  and a second one by ``secondary_domain`` which has to be declared
1716  first into the model.
1717  */
1718  size_type APIDECL add_linear_twodomain_term
1719  (model &md, const mesh_im &mim, const std::string &expr,
1720  size_type region, const std::string &secondary_domain,
1721  bool is_sym = false, bool is_coercive = false,
1722  const std::string &brickname = "", bool return_if_nonlin = false);
1723 
1724  /** Adds a nonlinear term given by a weak form language expression like
1725  ``add_nonlinear_term`` function but for an integration on a direct
1726  product of two domains, a first specfied by ``mim`` and ``region``
1727  and a second one by ``secondary_domain`` which has to be declared
1728  first into the model.
1729  */
1730  size_type APIDECL add_nonlinear_twodomain_term
1731  (model &md, const mesh_im &mim, const std::string &expr,
1732  size_type region, const std::string &secondary_domain,
1733  bool is_sym = false, bool is_coercive = false,
1734  const std::string &brickname = "");
1735 
1736  /** Adds a source term given by a weak form language expression like
1737  ``add_source_term`` function but for an integration on a direct
1738  product of two domains, a first specfied by ``mim`` and ``region``
1739  and a second one by ``secondary_domain`` which has to be declared
1740  first into the model.
1741  */
1742  size_type APIDECL add_twodomain_source_term
1743  (model &md, const mesh_im &mim, const std::string &expr,
1744  size_type region, const std::string &secondary_domain,
1745  const std::string &brickname = std::string(),
1746  const std::string &directvarname = std::string(),
1747  const std::string &directdataname = std::string(),
1748  bool return_if_nonlin = false);
1749 
1750 
1751  /** Add a Laplacian term on the variable `varname` (in fact with a minus :
1752  :math:`-\text{div}(\nabla u)`). If it is a vector
1753  valued variable, the Laplacian term is componentwise. `region` is an
1754  optional mesh region on which the term is added. Return the brick index
1755  in the model.
1756  */
1757  size_type APIDECL add_Laplacian_brick
1758  (model &md, const mesh_im &mim, const std::string &varname,
1759  size_type region = size_type(-1));
1760 
1761 
1762  /** Add an elliptic term on the variable `varname`.
1763  The shape of the elliptic
1764  term depends both on the variable and the data. This corresponds to a
1765  term $-\text{div}(a\nabla u)$ where $a$ is the data and $u$ the variable.
1766  The data can be a scalar, a matrix or an order four tensor. The variable
1767  can be vector valued or not. If the data is a scalar or a matrix and
1768  the variable is vector valued then the term is added componentwise.
1769  An order four tensor data is allowed for vector valued variable only.
1770  The data can be constant or describbed on a fem. Of course, when
1771  the data is a tensor describe on a finite element method (a tensor
1772  field) the data can be a huge vector. The components of the
1773  matrix/tensor have to be stored with the fortran order (columnwise) in
1774  the data vector (compatibility with blas). The symmetry and coercivity
1775  of the given matrix/tensor is not verified (but assumed). `region` is an
1776  optional mesh region on which the term is added. Note that for the real
1777  version which uses the high-level generic assembly language, `dataexpr`
1778  can be any regular expression of the high-level generic assembly
1779  language (like "1", "sin(X[0])" or "Norm(u)" for instance) even
1780  depending on model variables.
1781  Return the brick index in the model.
1782  */
1783  size_type APIDECL add_generic_elliptic_brick
1784  (model &md, const mesh_im &mim, const std::string &varname,
1785  const std::string &dataexpr, size_type region = size_type(-1));
1786 
1787 
1788  /** Add a source term on the variable `varname`. The source term is
1789  represented by `dataexpr` which could be any regular expression of the
1790  high-level generic assembly language (except for the complex version
1791  where it has to be a declared data of the model). `region` is an
1792  optional mesh region on which the term is
1793  added. An additional optional data `directdataname` can be provided. The
1794  corresponding data vector will be directly added to the right hand
1795  side without assembly. Return the brick index in the model.
1796  */
1797  size_type APIDECL add_source_term_brick
1798  (model &md, const mesh_im &mim, const std::string &varname,
1799  const std::string &dataexpr, size_type region = size_type(-1),
1800  const std::string &directdataname = std::string());
1801 
1802  /** Add a source term on the variable `varname` on a boundary `region`.
1803  The source term is
1804  represented by the data `dataepxpr` which could be any regular
1805  expression of the high-level generic assembly language (except
1806  for the complex version where it has to be a declared data of
1807  the model). A scalar product with the outward normal unit vector to
1808  the boundary is performed. The main aim of this brick is to represent
1809  a Neumann condition with a vector data without performing the
1810  scalar product with the normal as a pre-processing.
1811  */
1812  size_type APIDECL add_normal_source_term_brick
1813  (model &md, const mesh_im &mim, const std::string &varname,
1814  const std::string &dataexpr, size_type region);
1815 
1816 
1817  /** Add a (simple) Dirichlet condition on the variable `varname` and
1818  the mesh region `region`. The Dirichlet condition is prescribed by
1819  a simple post-treatment of the final linear system (tangent system
1820  for nonlinear problems) consisting of modifying the lines corresponding
1821  to the degree of freedom of the variable on `region` (0 outside the
1822  diagonal, 1 on the diagonal of the matrix and the expected value on
1823  the right hand side).
1824  The symmetry of the linear system is kept if all other bricks are
1825  symmetric.
1826  This brick is to be reserved for simple Dirichlet conditions (only dof
1827  declared on the correspodning boundary are prescribed). The application
1828  of this brick on reduced f.e.m. may be problematic. Intrinsic vectorial
1829  finite element method are not supported.
1830  `dataname` is the optional right hand side of the Dirichlet condition.
1831  It could be constant or (important) described on the same finite
1832  element method as `varname`.
1833  Returns the brick index in the model.
1834  */
1835  size_type APIDECL add_Dirichlet_condition_with_simplification
1836  (model &md, const std::string &varname, size_type region,
1837  const std::string &dataname = std::string());
1838 
1839 
1840  /** Add a Dirichlet condition on the variable `varname` and the mesh
1841  region `region`. This region should be a boundary. The Dirichlet
1842  condition is prescribed with a multiplier variable `multname` which
1843  should be first declared as a multiplier
1844  variable on the mesh region in the model. `dataname` is the optional
1845  right hand side of the Dirichlet condition. It could be constant or
1846  described on a fem; scalar or vector valued, depending on the variable
1847  on which the Dirichlet condition is prescribed. Return the brick index
1848  in the model.
1849  */
1850  size_type APIDECL add_Dirichlet_condition_with_multipliers
1851  (model &md, const mesh_im &mim, const std::string &varname,
1852  const std::string &multname, size_type region,
1853  const std::string &dataname = std::string());
1854 
1855  /** Same function as the previous one but the multipliers variable will
1856  be declared to the brick by the function. `mf_mult` is the finite element
1857  method on which the multiplier will be build (it will be restricted to
1858  the mesh region `region` and eventually some conflicting dofs with some
1859  other multiplier variables will be suppressed).
1860  */
1861  size_type APIDECL add_Dirichlet_condition_with_multipliers
1862  (model &md, const mesh_im &mim, const std::string &varname,
1863  const mesh_fem &mf_mult, size_type region,
1864  const std::string &dataname = std::string());
1865 
1866  /** Same function as the previous one but the `mf_mult` parameter is
1867  replaced by `degree`. The multiplier will be described on a standard
1868  finite element method of the corresponding degree.
1869  */
1870  size_type APIDECL add_Dirichlet_condition_with_multipliers
1871  (model &md, const mesh_im &mim, const std::string &varname,
1872  dim_type degree, size_type region,
1873  const std::string &dataname = std::string());
1874 
1875 
1876  /** When `ind_brick` is the index of a Dirichlet brick with multiplier on
1877  the model `md`, the function return the name of the multiplier variable.
1878  Otherwise, it has an undefined behavior.
1879  */
1880  const APIDECL std::string &mult_varname_Dirichlet(model &md, size_type ind_brick);
1881 
1882  /** Add a Dirichlet condition on the variable `varname` and the mesh
1883  region `region`. This region should be a boundary. The Dirichlet
1884  condition is prescribed with penalization. The penalization coefficient
1885  is intially `penalization_coeff` and will be added to the data of
1886  the model. `dataname` is the optional
1887  right hand side of the Dirichlet condition. It could be constant or
1888  described on a fem; scalar or vector valued, depending on the variable
1889  on which the Dirichlet condition is prescribed.
1890  `mf_mult` is an optional parameter which allows to weaken the
1891  Dirichlet condition specifying a multiplier space.
1892  Returns the brick index in the model.
1893  */
1894  size_type APIDECL add_Dirichlet_condition_with_penalization
1895  (model &md, const mesh_im &mim, const std::string &varname,
1896  scalar_type penalization_coeff, size_type region,
1897  const std::string &dataname = std::string(),
1898  const mesh_fem *mf_mult = 0);
1899 
1900  /** Add a Dirichlet condition on the variable `varname` and the mesh
1901  region `region`. This region should be a boundary. `Neumannterm`
1902  is the expression of the Neumann term (obtained by the Green formula)
1903  described as an expression of the high-level
1904  generic assembly language. This term can be obtained with
1905  md. Neumann_term(varname, region) once all volumic bricks have
1906  been added to the model. The Dirichlet
1907  condition is prescribed with Nitsche's method. `datag` is the optional
1908  right hand side of the Dirichlet condition. `datagamma0` is the
1909  Nitsche's method parameter. `theta` is a scalar value which can be
1910  positive or negative. `theta = 1` corresponds to the standard symmetric
1911  method which is conditionnaly coercive for `gamma0` small.
1912  `theta = -1` corresponds to the skew-symmetric method which is
1913  inconditionnaly coercive. `theta = 0` is the simplest method
1914  for which the second derivative of the Neumann term is not necessary
1915  even for nonlinear problems. Return the brick index in the model.
1916  */
1917  size_type APIDECL add_Dirichlet_condition_with_Nitsche_method
1918  (model &md, const mesh_im &mim, const std::string &varname,
1919  const std::string &Neumannterm,
1920  const std::string &datagamma0, size_type region,
1921  scalar_type theta = scalar_type(0),
1922  const std::string &datag = std::string());
1923 
1924 
1925  /** Add a Dirichlet condition to the normal component of the vector
1926  (or tensor) valued variable `varname` and the mesh
1927  region `region`. This region should be a boundary. The Dirichlet
1928  condition is prescribed with a multiplier variable `multname` which
1929  should be first declared as a multiplier
1930  variable on the mesh region in the model. `dataname` is the optional
1931  right hand side of the normal Dirichlet condition.
1932  It could be constant or
1933  described on a fem; scalar or vector valued, depending on the variable
1934  on which the Dirichlet condition is prescribed (scalar if the variable
1935  is vector valued, vector if the variable is tensor valued).
1936  Returns the brick index in the model.
1937  */
1938  size_type APIDECL add_normal_Dirichlet_condition_with_multipliers
1939  (model &md, const mesh_im &mim, const std::string &varname,
1940  const std::string &multname, size_type region,
1941  const std::string &dataname = std::string());
1942 
1943  /** Same function as the previous one but the multipliers variable will
1944  be declared to the brick by the function. `mf_mult` is the finite element
1945  method on which the multiplier will be build (it will be restricted to
1946  the mesh region `region` and eventually some conflicting dofs with some
1947  other multiplier variables will be suppressed).
1948  */
1949  size_type APIDECL add_normal_Dirichlet_condition_with_multipliers
1950  (model &md, const mesh_im &mim, const std::string &varname,
1951  const mesh_fem &mf_mult, size_type region,
1952  const std::string &dataname = std::string());
1953 
1954  /** Same function as the previous one but the `mf_mult` parameter is
1955  replaced by `degree`. The multiplier will be described on a standard
1956  finite element method of the corresponding degree.
1957  */
1958  size_type APIDECL add_normal_Dirichlet_condition_with_multipliers
1959  (model &md, const mesh_im &mim, const std::string &varname,
1960  dim_type degree, size_type region,
1961  const std::string &dataname = std::string());
1962 
1963  /** Add a Dirichlet condition to the normal component of the vector
1964  (or tensor) valued variable `varname` and the mesh
1965  region `region`. This region should be a boundary. The Dirichlet
1966  condition is prescribed with penalization. The penalization coefficient
1967  is intially `penalization_coeff` and will be added to the data of
1968  the model. `dataname` is the optional
1969  right hand side of the Dirichlet condition. It could be constant or
1970  described on a fem; scalar or vector valued, depending on the variable
1971  on which the Dirichlet condition is prescribed (scalar if the variable
1972  is vector valued, vector if the variable is tensor valued).
1973  `mf_mult` is an optional parameter which allows to weaken the
1974  Dirichlet condition specifying a multiplier space.
1975  Return the brick index in the model.
1976  */
1977  size_type APIDECL add_normal_Dirichlet_condition_with_penalization
1978  (model &md, const mesh_im &mim, const std::string &varname,
1979  scalar_type penalization_coeff, size_type region,
1980  const std::string &dataname = std::string(),
1981  const mesh_fem *mf_mult = 0);
1982 
1983 
1984 
1985  /** Add a Dirichlet condition on the normal component of the variable
1986  `varname` and the mesh
1987  region `region`. This region should be a boundary. `Neumannterm`
1988  is the expression of the Neumann term (obtained by the Green formula)
1989  described as an expression of the high-level
1990  generic assembly language. This term can be obtained with
1991  md.Neumann_term(varname, region) once all volumic bricks have
1992  been added to the model.The Dirichlet
1993  condition is prescribed with Nitsche's method. `datag` is the optional
1994  scalar right hand side of the Dirichlet condition. `datagamma0` is the
1995  Nitsche's method parameter. `theta` is a scalar value which can be
1996  positive or negative. `theta = 1` corresponds to the standard symmetric
1997  method which is conditionnaly coercive for `gamma0` small.
1998  `theta = -1` corresponds to the skew-symmetric method which is
1999  inconditionnaly coercive. `theta = 0` is the simplest method
2000  for which the second derivative of the Neumann term is not necessary
2001  even for nonlinear problems. Return the brick index in the model.
2002  */
2003  size_type APIDECL add_normal_Dirichlet_condition_with_Nitsche_method
2004  (model &md, const mesh_im &mim, const std::string &varname,
2005  const std::string &Neumannterm, const std::string &datagamma0,
2006  size_type region, scalar_type theta = scalar_type(0),
2007  const std::string &datag = std::string());
2008 
2009 
2010  /** Add some pointwise constraints on the variable `varname` thanks to
2011  a penalization. The penalization coefficient is initially
2012  `penalization_coeff` and will be added to the data of the model.
2013  The conditions are prescribed on a set of points given in the data
2014  `dataname_pt` whose dimension is the number of points times the dimension
2015  of the mesh. If the variable represents a vector field, the data
2016  `dataname_unitv` represents a vector of dimension the number of points
2017  times the dimension of the vector field which should store some
2018  unit vectors. In that case the prescribed constraint is the scalar
2019  product of the variable at the corresponding point with the corresponding
2020  unit vector.
2021  The optional data `dataname_val` is the vector of values to be prescribed
2022  at the different points.
2023  This brick is specifically designed to kill rigid displacement
2024  in a Neumann problem.
2025  */
2026  size_type APIDECL add_pointwise_constraints_with_penalization
2027  (model &md, const std::string &varname,
2028  scalar_type penalisation_coeff, const std::string &dataname_pt,
2029  const std::string &dataname_unitv = std::string(),
2030  const std::string &dataname_val = std::string());
2031 
2032 
2033  /** Add some pointwise constraints on the variable `varname` using a given
2034  multiplier `multname`.
2035  The conditions are prescribed on a set of points given in the data
2036  `dataname_pt` whose dimension is the number of points times the dimension
2037  of the mesh.
2038  The multiplier variable should be a fixed size variable of size the
2039  number of points.
2040  If the variable represents a vector field, the data
2041  `dataname_unitv` represents a vector of dimension the number of points
2042  times the dimension of the vector field which should store some
2043  unit vectors. In that case the prescribed constraint is the scalar
2044  product of the variable at the corresponding point with the corresponding
2045  unit vector.
2046  The optional data `dataname_val` is the vector of values to be prescribed
2047  at the different points.
2048  This brick is specifically designed to kill rigid displacement
2049  in a Neumann problem.
2050  */
2051  size_type APIDECL add_pointwise_constraints_with_given_multipliers
2052  (model &md, const std::string &varname,
2053  const std::string &multname, const std::string &dataname_pt,
2054  const std::string &dataname_unitv = std::string(),
2055  const std::string &dataname_val = std::string());
2056 
2057  /** Add some pointwise constraints on the variable `varname` using
2058  multiplier. The multiplier variable is automatically added to the model.
2059  The conditions are prescribed on a set of points given in the data
2060  `dataname_pt` whose dimension is the number of points times the dimension
2061  of the mesh.
2062  If the variable represents a vector field, the data
2063  `dataname_unitv` represents a vector of dimension the number of points
2064  times the dimension of the vector field which should store some
2065  unit vectors. In that case the prescribed constraint is the scalar
2066  product of the variable at the corresponding point with the corresponding
2067  unit vector.
2068  The optional data `dataname_val` is the vector of values to be prescribed
2069  at the different points.
2070  This brick is specifically designed to kill rigid displacement
2071  in a Neumann problem.
2072  */
2073  size_type APIDECL add_pointwise_constraints_with_multipliers
2074  (model &md, const std::string &varname, const std::string &dataname_pt,
2075  const std::string &dataname_unitv = std::string(),
2076  const std::string &dataname_val = std::string());
2077 
2078 
2079  /** Change the penalization coefficient of a Dirichlet condition with
2080  penalization brick. If the brick is not of this kind,
2081  this function has an undefined behavior.
2082  */
2083  void APIDECL change_penalization_coeff(model &md, size_type ind_brick,
2084  scalar_type penalisation_coeff);
2085 
2086  /** Add a generalized Dirichlet condition on the variable `varname` and
2087  the mesh region `region`. This version is for vector field.
2088  It prescribes a condition @f$ Hu = r @f$ where `H` is a matrix field.
2089  This region should be a boundary. The Dirichlet
2090  condition is prescribed with a multiplier variable `multname` which
2091  should be first declared as a multiplier
2092  variable on the mesh region in the model. `dataname` is the
2093  right hand side of the Dirichlet condition. It could be constant or
2094  described on a fem; scalar or vector valued, depending on the variable
2095  on which the Dirichlet condition is prescribed. `Hname' is the data
2096  corresponding to the matrix field `H`. It has to be a constant matrix
2097  or described on a scalar fem. Return the brick index in the model.
2098  */
2099  size_type APIDECL add_generalized_Dirichlet_condition_with_multipliers
2100  (model &md, const mesh_im &mim, const std::string &varname,
2101  const std::string &multname, size_type region,
2102  const std::string &dataname, const std::string &Hname);
2103 
2104  /** Same function as the preceeding one but the multipliers variable will
2105  be declared to the brick by the function. `mf_mult` is the finite element
2106  method on which the multiplier will be build (it will be restricted to
2107  the mesh region `region` and eventually some conflicting dofs with some
2108  other multiplier variables will be suppressed).
2109  */
2110  size_type APIDECL add_generalized_Dirichlet_condition_with_multipliers
2111  (model &md, const mesh_im &mim, const std::string &varname,
2112  const mesh_fem &mf_mult, size_type region,
2113  const std::string &dataname, const std::string &Hname);
2114 
2115  /** Same function as the preceeding one but the `mf_mult` parameter is
2116  replaced by `degree`. The multiplier will be described on a standard
2117  finite element method of the corresponding degree.
2118  */
2119  size_type APIDECL add_generalized_Dirichlet_condition_with_multipliers
2120  (model &md, const mesh_im &mim, const std::string &varname,
2121  dim_type degree, size_type region,
2122  const std::string &dataname, const std::string &Hname);
2123 
2124  /** Add a Dirichlet condition on the variable `varname` and the mesh
2125  region `region`. This version is for vector field.
2126  It prescribes a condition @f$ Hu = r @f$ where `H` is a matrix field.
2127  This region should be a boundary. This region should be a boundary.
2128  The Dirichlet
2129  condition is prescribed with penalization. The penalization coefficient
2130  is intially `penalization_coeff` and will be added to the data of
2131  the model. `dataname` is the
2132  right hand side of the Dirichlet condition. It could be constant or
2133  described on a fem; scalar or vector valued, depending on the variable
2134  on which the Dirichlet condition is prescribed. `Hname' is the data
2135  corresponding to the matrix field `H`. It has to be a constant matrix
2136  or described on a scalar fem. `mf_mult` is an optional parameter
2137  which allows to weaken the Dirichlet condition specifying a
2138  multiplier space. Return the brick index in the model.
2139  */
2140  size_type APIDECL add_generalized_Dirichlet_condition_with_penalization
2141  (model &md, const mesh_im &mim, const std::string &varname,
2142  scalar_type penalization_coeff, size_type region,
2143  const std::string &dataname, const std::string &Hname,
2144  const mesh_fem *mf_mult = 0);
2145 
2146  /** Add a Dirichlet condition on the variable `varname` and the mesh
2147  region `region`. This region should be a boundary. This version
2148  is for vector field. It prescribes a condition
2149  @f$ Hu = r @f$ where `H` is a matrix field. `Neumannterm`
2150  is the expression of the Neumann term (obtained by the Green formula)
2151  described as an expression of the high-level
2152  generic assembly language. of the high-level
2153  generic assembly language. This term can be obtained with
2154  md.Neumann_term(varname, region) once all volumic bricks have
2155  been added to the model. The Dirichlet
2156  condition is prescribed with Nitsche's method. `datag` is the optional
2157  right hand side of the Dirichlet condition. `datagamma0` is the
2158  Nitsche's method parameter. `theta` is a scalar value which can be
2159  positive or negative. `theta = 1` corresponds to the standard symmetric
2160  method which is conditionnaly coercive for `gamma0` small.
2161  `theta = -1` corresponds to the skew-symmetric method which is
2162  inconditionnaly coercive. `theta = 0` is the simplest method
2163  for which the second derivative of the Neumann term is not necessary
2164  even for nonlinear problems. Return the brick index in the model.
2165  */
2166  size_type APIDECL add_generalized_Dirichlet_condition_with_Nitsche_method
2167  (model &md, const mesh_im &mim, const std::string &varname,
2168  const std::string &Neumannterm, const std::string &datagamma0,
2169  size_type region, scalar_type theta,
2170  const std::string &datag, const std::string &dataH);
2171 
2172 
2173  /** Add a Helmoltz brick to the model. This corresponds to the scalar
2174  equation (@f$\Delta u + k^2u = 0@f$, with @f$K=k^2@f$).
2175  The weak formulation is (@f$\int k^2 u.v - \nabla u.\nabla v@f$)
2176 
2177  `dataexpr` should contain the wave number $k$. It can be real or
2178  complex.
2179  */
2180  size_type APIDECL add_Helmholtz_brick(model &md, const mesh_im &mim,
2181  const std::string &varname,
2182  const std::string &dataexpr,
2183  size_type region = size_type(-1));
2184 
2185 
2186  /** Add a Fourier-Robin brick to the model. This correspond to the weak term
2187  (@f$\int (qu).v @f$) on a boundary. It is used to represent a
2188  Fourier-Robin boundary condition.
2189 
2190  `dataexpr` is the parameter $q$ which should be a
2191  (@f$N\times N@f$) matrix term, where $N$ is the dimension of the
2192  variable `varname`. It can be an arbitrary valid expression of the
2193  high-level generic assembly language (except for the complex version
2194  for which it should be a data of the model). Note that an additional
2195  right hand side can be added with a source term brick.
2196  */
2197  size_type APIDECL add_Fourier_Robin_brick(model &md, const mesh_im &mim,
2198  const std::string &varname,
2199  const std::string &dataexpr,
2200  size_type region);
2201 
2202 
2203  // Constraint brick.
2204  model_real_sparse_matrix APIDECL &set_private_data_brick_real_matrix
2205  (model &md, size_type indbrick);
2206  model_real_plain_vector APIDECL &set_private_data_brick_real_rhs
2207  (model &md, size_type indbrick);
2208  model_complex_sparse_matrix APIDECL &set_private_data_brick_complex_matrix
2209  (model &md, size_type indbrick);
2210  model_complex_plain_vector APIDECL &set_private_data_brick_complex_rhs
2211  (model &md, size_type indbrick);
2212  size_type APIDECL add_constraint_with_penalization
2213  (model &md, const std::string &varname, scalar_type penalisation_coeff);
2214  size_type APIDECL add_constraint_with_multipliers
2215  (model &md, const std::string &varname, const std::string &multname);
2216 
2217  void set_private_data_rhs
2218  (model &md, size_type indbrick, const std::string &varname);
2219 
2220  template <typename VECT, typename T>
2221  void set_private_data_rhs(model &md, size_type ind,
2222  const VECT &L, T) {
2223  model_real_plain_vector &LL = set_private_data_brick_real_rhs(md, ind);
2224  gmm::resize(LL, gmm::vect_size(L));
2225  gmm::copy(L, LL);
2226  }
2227 
2228  template <typename VECT, typename T>
2229  void set_private_data_rhs(model &md, size_type ind, const VECT &L,
2230  std::complex<T>) {
2231  model_complex_plain_vector &LL = set_private_data_brick_complex_rhs(md, ind);
2232  gmm::resize(LL, gmm::vect_size(L));
2233  gmm::copy(L, LL);
2234  }
2235 
2236  /** For some specific bricks having an internal right hand side vector
2237  (explicit bricks: 'constraint brick' and 'explicit rhs brick'),
2238  set this rhs.
2239  */
2240  template <typename VECT>
2241  void set_private_data_rhs(model &md, size_type indbrick, const VECT &L) {
2242  typedef typename gmm::linalg_traits<VECT>::value_type T;
2243  set_private_data_rhs(md, indbrick, L, T());
2244  }
2245 
2246  template <typename MAT, typename T>
2247  void set_private_data_matrix(model &md, size_type ind,
2248  const MAT &B, T) {
2249  model_real_sparse_matrix &BB = set_private_data_brick_real_matrix(md, ind);
2250  gmm::resize(BB, gmm::mat_nrows(B), gmm::mat_ncols(B));
2251  gmm::copy(B, BB);
2252  }
2253 
2254  template <typename MAT, typename T>
2255  void set_private_data_matrix(model &md, size_type ind, const MAT &B,
2256  std::complex<T>) {
2257  model_complex_sparse_matrix &BB
2258  = set_private_data_brick_complex_matrix(md, ind);
2259  gmm::resize(BB, gmm::mat_nrows(B), gmm::mat_ncols(B));
2260  gmm::copy(B, BB);
2261  }
2262 
2263  /** For some specific bricks having an internal sparse matrix
2264  (explicit bricks: 'constraint brick' and 'explicit matrix brick'),
2265  set this matrix. @*/
2266  template <typename MAT>
2267  void set_private_data_matrix(model &md, size_type indbrick,
2268  const MAT &B) {
2269  typedef typename gmm::linalg_traits<MAT>::value_type T;
2270  set_private_data_matrix(md, indbrick, B, T());
2271  }
2272 
2273  /** Add an additional explicit penalized constraint on the variable
2274  `varname`. The constraint is $BU=L$ with `B` being a rectangular
2275  sparse matrix.
2276  Be aware that `B` should not contain a plain row, otherwise the whole
2277  tangent matrix will be plain. It is possible to change the constraint
2278  at any time with the methods set_private_matrix and set_private_rhs.
2279  The method change_penalization_coeff can also be used.
2280  */
2281  template <typename MAT, typename VECT>
2282  size_type add_constraint_with_penalization
2283  (model &md, const std::string &varname, scalar_type penalisation_coeff,
2284  const MAT &B, const VECT &L) {
2285  size_type ind
2286  = add_constraint_with_penalization(md, varname, penalisation_coeff);
2287  size_type n = gmm::mat_nrows(B), m = gmm::mat_ncols(B);
2288  set_private_data_rhs(md, ind, L);
2289  set_private_data_matrix(md, ind, B);
2290  return ind;
2291  }
2292 
2293  /** Add an additional explicit constraint on the variable `varname` thank to
2294  a multiplier `multname` peviously added to the model (should be a fixed
2295  size variable).
2296  The constraint is $BU=L$ with `B` being a rectangular sparse matrix.
2297  It is possible to change the constraint
2298  at any time with the methods set_private_matrix
2299  and set_private_rhs.
2300  */
2301  template <typename MAT, typename VECT>
2302  size_type add_constraint_with_multipliers
2303  (model &md, const std::string &varname, const std::string &multname,
2304  const MAT &B, const VECT &L) {
2305  size_type ind = add_constraint_with_multipliers(md, varname, multname);
2306  set_private_data_rhs(md, ind, L);
2307  set_private_data_matrix(md, ind, B);
2308  return ind;
2309  }
2310 
2311  template <typename MAT>
2312  size_type add_constraint_with_multipliers
2313  (model &md, const std::string &varname, const std::string &multname,
2314  const MAT &B, const std::string &Lname) {
2315  size_type ind = add_constraint_with_multipliers(md, varname, multname);
2316  set_private_data_rhs(md, ind, Lname);
2317  set_private_data_matrix(md, ind, B);
2318  return ind;
2319  }
2320 
2321  size_type APIDECL add_explicit_matrix(model &md, const std::string &varname1,
2322  const std::string &varname2,
2323  bool issymmetric, bool iscoercive);
2324  size_type APIDECL add_explicit_rhs(model &md, const std::string &varname);
2325 
2326  /** Add a brick reprenting an explicit matrix to be added to the tangent
2327  linear system relatively to the variables 'varname1' and 'varname2'.
2328  The given matrix should have as many rows as the dimension of
2329  'varname1' and as many columns as the dimension of 'varname2'.
2330  If the two variables are different and if `issymmetric' is set to true
2331  then the transpose of the matrix is also added to the tangent system
2332  (default is false). set `iscoercive` to true if the term does not
2333  affect the coercivity of the tangent system (default is false).
2334  The matrix can be changed by the command set_private_matrix.
2335  */
2336  template <typename MAT>
2337  size_type add_explicit_matrix(model &md, const std::string &varname1,
2338  const std::string &varname2, const MAT &B,
2339  bool issymmetric = false,
2340  bool iscoercive = false) {
2341  size_type ind = add_explicit_matrix(md, varname1, varname2,
2342  issymmetric, iscoercive);
2343  set_private_data_matrix(md, ind, B);
2344  return ind;
2345  }
2346 
2347  /** Add a brick representing an explicit right hand side to be added to
2348  the right hand side of the tangent
2349  linear system relatively to the variable 'varname'.
2350  The given rhs should have the same size than the dimension of
2351  'varname'. The rhs can be changed by the command set_private_rhs.
2352  */
2353  template <typename VECT>
2354  size_type add_explicit_rhs(model &md, const std::string &varname,
2355  const VECT &L) {
2356  size_type ind = add_explicit_rhs(md, varname);
2357  set_private_data_rhs(md, ind, L);
2358  return ind;
2359  }
2360 
2361 
2362  /** Linear elasticity brick ( @f$ \int \sigma(u):\varepsilon(v) @f$ ).
2363  for isotropic material. Parametrized by the Lamé coefficients
2364  lambda and mu.
2365  */
2366  size_type APIDECL add_isotropic_linearized_elasticity_brick
2367  (model &md, const mesh_im &mim, const std::string &varname,
2368  const std::string &dataname_lambda, const std::string &dataname_mu,
2369  size_type region = size_type(-1),
2370  const std::string &dataname_preconstraint = std::string());
2371 
2372  /** Linear elasticity brick ( @f$ \int \sigma(u):\varepsilon(v) @f$ ).
2373  for isotropic material. Parametrized by Young modulus and Poisson ratio
2374  For two-dimensional problems, corresponds to the plane strain
2375  approximation
2376  ( @f$ \lambda = E\nu/((1+\nu)(1-2\nu)), \mu = E/(2(1+\nu)) @f$ ).
2377  Corresponds to the standard model for three-dimensional problems.
2378  */
2379  size_type APIDECL add_isotropic_linearized_elasticity_pstrain_brick
2380  (model &md, const mesh_im &mim, const std::string &varname,
2381  const std::string &data_E, const std::string &data_nu,
2382  size_type region);
2383 
2384  /**
2385  Linear elasticity brick ( @f$ \int \sigma(u):\varepsilon(v) @f$ ).
2386  for isotropic material. Parametrized by Young modulus and Poisson ratio.
2387  For two-dimensional problems, corresponds to the plane stress
2388  approximation
2389  ( @f$ \lambda^* = E\nu/(1-\nu^2), \mu = E/(2(1+\nu)) @f$ ).
2390  Corresponds to the standard model for three-dimensional problems.
2391  */
2392  size_type APIDECL add_isotropic_linearized_elasticity_pstress_brick
2393  (model &md, const mesh_im &mim, const std::string &varname,
2394  const std::string &data_E, const std::string &data_nu,
2395  size_type region);
2396 
2397  void APIDECL compute_isotropic_linearized_Von_Mises_or_Tresca
2398  (model &md, const std::string &varname, const std::string &dataname_lambda,
2399  const std::string &dataname_mu, const mesh_fem &mf_vm,
2400  model_real_plain_vector &VM, bool tresca);
2401 
2402  /**
2403  Compute the Von-Mises stress or the Tresca stress of a field
2404  (only valid for isotropic linearized elasticity in 3D)
2405  Parametrized by Lame coefficients.
2406  */
2407  template <class VECTVM>
2408  void compute_isotropic_linearized_Von_Mises_or_Tresca
2409  (model &md, const std::string &varname, const std::string &dataname_lambda,
2410  const std::string &dataname_mu, const mesh_fem &mf_vm,
2411  VECTVM &VM, bool tresca) {
2412  model_real_plain_vector VMM(mf_vm.nb_dof());
2413  compute_isotropic_linearized_Von_Mises_or_Tresca
2414  (md, varname, dataname_lambda, dataname_mu, mf_vm, VMM, tresca);
2415  gmm::copy(VMM, VM);
2416  }
2417 
2418  /**
2419  Compute the Von-Mises stress of a displacement field for isotropic
2420  linearized elasticity in 3D or in 2D with plane strain assumption.
2421  Parametrized by Young modulus and Poisson ratio.
2422  */
2423  void APIDECL compute_isotropic_linearized_Von_Mises_pstrain
2424  (model &md, const std::string &varname, const std::string &data_E,
2425  const std::string &data_nu, const mesh_fem &mf_vm,
2426  model_real_plain_vector &VM);
2427 
2428  /**
2429  Compute the Von-Mises stress of a displacement field for isotropic
2430  linearized elasticity in 3D or in 2D with plane stress assumption.
2431  Parametrized by Young modulus and Poisson ratio.
2432  */
2433  void APIDECL compute_isotropic_linearized_Von_Mises_pstress
2434  (model &md, const std::string &varname, const std::string &data_E,
2435  const std::string &data_nu, const mesh_fem &mf_vm,
2436  model_real_plain_vector &VM);
2437 
2438 
2439  /**
2440  Mixed linear incompressibility condition brick.
2441 
2442  Update the tangent matrix with a pressure term:
2443  @f[
2444  T \longrightarrow
2445  \begin{array}{ll} T & B \\ B^t & M \end{array}
2446  @f]
2447  with @f$ B = - \int p.div u @f$ and
2448  @f$ M = \int \epsilon p.q @f$ ( @f$ \epsilon @f$ is an optional
2449  penalization coefficient).
2450 
2451  Be aware that an inf-sup condition between the finite element method
2452  describing the rpressure and the primal variable has to be satisfied.
2453 
2454  For nearly incompressible elasticity,
2455  @f[ p = -\lambda \textrm{div}~u @f]
2456  @f[ \sigma = 2 \mu \varepsilon(u) -p I @f]
2457  */
2458  size_type APIDECL add_linear_incompressibility
2459  (model &md, const mesh_im &mim, const std::string &varname,
2460  const std::string &multname_pressure, size_type region = size_type(-1),
2461  const std::string &dataexpr_penal_coeff = std::string());
2462 
2463  /** Mass brick ( @f$ \int \rho u.v @f$ ).
2464  Add a mass matix on a variable (eventually with a specified region).
2465  If the parameter $\rho$ is omitted it is assumed to be equal to 1.
2466  */
2467  size_type APIDECL add_mass_brick
2468  (model &md, const mesh_im &mim, const std::string &varname,
2469  const std::string &dataexpr_rho = std::string(),
2470  size_type region = size_type(-1));
2471 
2472  /** Lumped mass brick for first order.
2473  Add a lumped mass matix for first order on a variable (eventually with a specified region).
2474  If the parameter $\rho$ is omitted it is assumed to be equal to 1.
2475  */
2476  size_type APIDECL add_lumped_mass_for_first_order_brick
2477  (model &md, const mesh_im &mim, const std::string &varname,
2478  const std::string &dataexpr_rho = std::string(),
2479  size_type region = size_type(-1));
2480 
2481  /** Basic d/dt brick ( @f$ \int \rho ((u^{n+1}-u^n)/dt).v @f$ ).
2482  Add the standard discretization of a first order time derivative. The
2483  parameter $rho$ is the density which could be omitted (the defaul value
2484  is 1). This brick should be used in addition to a time dispatcher for the
2485  other terms.
2486  */
2487  size_type APIDECL add_basic_d_on_dt_brick
2488  (model &md, const mesh_im &mim, const std::string &varname,
2489  const std::string &dataname_dt,
2490  const std::string &dataname_rho = std::string(),
2491  size_type region = size_type(-1));
2492 
2493  /** Basic d2/dt2 brick ( @f$ \int \rho ((u^{n+1}-u^n)/(\alpha dt^2) - v^n/(\alpha dt) ).w @f$ ).
2494  Add the standard discretization of a second order time derivative. The
2495  parameter $rho$ is the density which could be omitted (the defaul value
2496  is 1). This brick should be used in addition to a time dispatcher for the
2497  other terms. The time derivative $v$ of the variable $u$ is preferably
2498  computed as a post-traitement which depends on each scheme.
2499  */
2500  size_type APIDECL add_basic_d2_on_dt2_brick
2501  (model &md, const mesh_im &mim, const std::string &varnameU,
2502  const std::string &datanameV,
2503  const std::string &dataname_dt,
2504  const std::string &dataname_alpha,
2505  const std::string &dataname_rho = std::string(),
2506  size_type region = size_type(-1));
2507 
2508 
2509 } /* end of namespace getfem. */
2510 
2511 
2512 #endif /* GETFEM_MODELS_H_*/
getfem::model::is_complex
bool is_complex() const
Boolean which says if the model deals with real or complex unknowns and data.
Definition: getfem_models.h:565
getfem::change_penalization_coeff
void APIDECL change_penalization_coeff(model &md, size_type ind_brick, scalar_type penalisation_coeff)
Change the penalization coefficient of a Dirichlet condition with penalization brick.
Definition: getfem_models.cc:4842
getfem::model::add_elementary_transformation
void add_elementary_transformation(const std::string &name, pelementary_transformation ptrans)
Add an elementary transformation to the model to be used with the generic assembly.
Definition: getfem_models.h:1135
getfem::add_nonlinear_twodomain_term
size_type APIDECL add_nonlinear_twodomain_term(model &md, const mesh_im &mim, const std::string &expr, size_type region, const std::string &secondary_domain, bool is_sym=false, bool is_coercive=false, const std::string &brickname="")
Adds a nonlinear term given by a weak form language expression like add_nonlinear_term function but f...
Definition: getfem_models.cc:3685
getfem::mult_varname_Dirichlet
const APIDECL std::string & mult_varname_Dirichlet(model &md, size_type ind_brick)
When ind_brick is the index of a Dirichlet brick with multiplier on the model md, the function return...
Definition: getfem_models.cc:4705
getfem::virtual_dispatcher
The time dispatcher object modify the result of a brick in order to apply a time integration scheme.
Definition: getfem_models.h:1263
getfem::add_normal_Dirichlet_condition_with_multipliers
size_type APIDECL add_normal_Dirichlet_condition_with_multipliers(model &md, const mesh_im &mim, const std::string &varname, const std::string &multname, size_type region, const std::string &dataname=std::string())
Add a Dirichlet condition to the normal component of the vector (or tensor) valued variable varname a...
Definition: getfem_models.cc:4730
getfem::compute_isotropic_linearized_Von_Mises_pstrain
void APIDECL compute_isotropic_linearized_Von_Mises_pstrain(model &md, const std::string &varname, const std::string &data_E, const std::string &data_nu, const mesh_fem &mf_vm, model_real_plain_vector &VM)
Compute the Von-Mises stress of a displacement field for isotropic linearized elasticity in 3D or in ...
Definition: getfem_models.cc:6256
getfem::add_generalized_Dirichlet_condition_with_penalization
size_type APIDECL add_generalized_Dirichlet_condition_with_penalization(model &md, const mesh_im &mim, const std::string &varname, scalar_type penalization_coeff, size_type region, const std::string &dataname, const std::string &Hname, const mesh_fem *mf_mult=0)
Add a Dirichlet condition on the variable varname and the mesh region region.
Definition: getfem_models.cc:4821
getfem::add_linear_twodomain_term
size_type APIDECL add_linear_twodomain_term(model &md, const mesh_im &mim, const std::string &expr, size_type region, const std::string &secondary_domain, bool is_sym=false, bool is_coercive=false, const std::string &brickname="", bool return_if_nonlin=false)
Adds a linear term given by a weak form language expression like add_linear_term function but for an ...
Definition: getfem_models.cc:3594
getfem::model::macro_exists
bool macro_exists(const std::string &name) const
Says if a macro of that name has been defined.
Definition: getfem_models.h:901
getfem::model::elementary_transformation
pelementary_transformation elementary_transformation(const std::string &name) const
Get a pointer to the elementary transformation name.
Definition: getfem_models.h:1143
getfem::add_midpoint_dispatcher
void APIDECL add_midpoint_dispatcher(model &md, dal::bit_vector ibricks)
Add a midpoint time dispatcher to a list of bricks.
Definition: getfem_models.cc:7350
gmm::resize
void resize(M &v, size_type m, size_type n)
*‍/
Definition: gmm_blas.h:231
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
getfem::model::is_linear
bool is_linear() const
Return true if all the model terms are linear.
Definition: getfem_models.h:583
getfem::model::disable_brick
void disable_brick(size_type ib)
Disable a brick.
Definition: getfem_models.h:515
getfem::mesh_im
Describe an integration method linked to a mesh.
Definition: getfem_mesh_im.h:47
getfem::add_pointwise_constraints_with_penalization
size_type APIDECL add_pointwise_constraints_with_penalization(model &md, const std::string &varname, scalar_type penalisation_coeff, const std::string &dataname_pt, const std::string &dataname_unitv=std::string(), const std::string &dataname_val=std::string())
Add some pointwise constraints on the variable varname thanks to a penalization.
Definition: getfem_models.cc:5362
getfem::add_isotropic_linearized_elasticity_pstrain_brick
size_type APIDECL add_isotropic_linearized_elasticity_pstrain_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &data_E, const std::string &data_nu, size_type region)
Linear elasticity brick ( ).
Definition: getfem_models.cc:6142
getfem::add_Dirichlet_condition_with_multipliers
size_type APIDECL add_Dirichlet_condition_with_multipliers(model &md, const mesh_im &mim, const std::string &varname, const std::string &multname, size_type region, const std::string &dataname=std::string())
Add a Dirichlet condition on the variable varname and the mesh region region.
Definition: getfem_models.cc:4671
getfem::add_linear_term
size_type APIDECL add_linear_term(model &md, const mesh_im &mim, const std::string &expr, size_type region=size_type(-1), bool is_sym=false, bool is_coercive=false, const std::string &brickname="", bool return_if_nonlin=false)
Add a term given by the weak form language expression expr which will be assembled in region region a...
Definition: getfem_models.cc:3586
getfem::model::add_initialized_fem_data
void add_initialized_fem_data(const std::string &name, const mesh_fem &mf, const VECT &v)
Add an initialized fixed size data to the model, assumed to be a vector field if the size of the vect...
Definition: getfem_models.h:835
getfem::model::nb_primary_dof
size_type nb_primary_dof() const
Number of primary degrees of freedom in the model.
Definition: getfem_models.h:592
getfem::add_nonlinear_term
size_type APIDECL add_nonlinear_term(model &md, const mesh_im &mim, const std::string &expr, size_type region=size_type(-1), bool is_sym=false, bool is_coercive=false, const std::string &brickname="")
Add a nonlinear term given by the weak form language expression expr which will be assembled in regio...
Definition: getfem_models.cc:3678
getfem::model::nb_internal_dof
size_type nb_internal_dof() const
Number of internal degrees of freedom in the model.
Definition: getfem_models.h:589
getfem::add_Dirichlet_condition_with_Nitsche_method
size_type APIDECL add_Dirichlet_condition_with_Nitsche_method(model &md, const mesh_im &mim, const std::string &varname, const std::string &Neumannterm, const std::string &datagamma0, size_type region, scalar_type theta=scalar_type(0), const std::string &datag=std::string())
Add a Dirichlet condition on the variable varname and the mesh region region.
Definition: getfem_models.cc:5040
getfem::model::enable_brick
void enable_brick(size_type ib)
Enable a brick.
Definition: getfem_models.h:521
getfem::add_generic_elliptic_brick
size_type APIDECL add_generic_elliptic_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataexpr, size_type region=size_type(-1))
Add an elliptic term on the variable varname.
Definition: getfem_models.cc:3946
getfem::add_theta_method_dispatcher
void APIDECL add_theta_method_dispatcher(model &md, dal::bit_vector ibricks, const std::string &THETA)
Add a theta-method time dispatcher to a list of bricks.
Definition: getfem_models.cc:7044
getfem::add_generalized_Dirichlet_condition_with_multipliers
size_type APIDECL add_generalized_Dirichlet_condition_with_multipliers(model &md, const mesh_im &mim, const std::string &varname, const std::string &multname, size_type region, const std::string &dataname, const std::string &Hname)
Add a generalized Dirichlet condition on the variable varname and the mesh region region.
Definition: getfem_models.cc:4785
getfem::add_mass_brick
size_type APIDECL add_mass_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataexpr_rho=std::string(), size_type region=size_type(-1))
Mass brick ( ).
Definition: getfem_models.cc:6526
getfem::add_linear_incompressibility
size_type APIDECL add_linear_incompressibility(model &md, const mesh_im &mim, const std::string &varname, const std::string &multname_pressure, size_type region=size_type(-1), const std::string &dataexpr_penal_coeff=std::string())
Mixed linear incompressibility condition brick.
Definition: getfem_models.cc:6377
getfem::velocity_update_for_order_two_theta_method
void APIDECL velocity_update_for_order_two_theta_method(model &md, const std::string &U, const std::string &V, const std::string &pdt, const std::string &ptheta)
Function which udpate the velocity $v^{n+1}$ after the computation of the displacement $u^{n+1}$ and ...
Definition: getfem_models.cc:7051
getfem::model::add_initialized_fixed_size_data
void add_initialized_fixed_size_data(const std::string &name, const VECT &v, const bgeot::multi_index &sizes)
Add a fixed size data (assumed to be a vector) to the model and initialized with v.
Definition: getfem_models.h:755
getfem::mesh_fem
Describe a finite element method linked to a mesh.
Definition: getfem_mesh_fem.h:148
getfem_partial_mesh_fem.h
a subclass of getfem::mesh_fem which allows to eliminate a number of dof of the original mesh_fem.
getfem::PREFIX_OLD
const auto PREFIX_OLD
A prefix to refer to the previous version of a variable.
Definition: getfem_models.h:98
getfem::model::add_initialized_fem_data
void add_initialized_fem_data(const std::string &name, const mesh_fem &mf, const VECT &v, const bgeot::multi_index &sizes)
Add a fixed size data to the model.
Definition: getfem_models.h:848
getfem::add_Helmholtz_brick
size_type APIDECL add_Helmholtz_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataexpr, size_type region=size_type(-1))
Add a Helmoltz brick to the model.
Definition: getfem_models.cc:5516
getfem::model
`‘Model’' variables store the variables, the data and the description of a model.
Definition: getfem_models.h:114
getfem
GEneric Tool for Finite Element Methods.
Definition: getfem_accumulated_distro.h:46
getfem::model::complex_rhs
const model_complex_plain_vector & complex_rhs() const
Gives access to the right hand side of the tangent linear system.
Definition: getfem_models.h:981
getfem::virtual_brick::complex_post_assembly_in_serial
virtual void complex_post_assembly_in_serial(const model &, size_type, const model::varnamelist &, const model::varnamelist &, const model::mimlist &, model::complex_matlist &, model::complex_veclist &, model::complex_veclist &, size_type, build_version) const
Peform any post assembly action for complex terms.
Definition: getfem_models.h:1579
getfem::virtual_brick::complex_pre_assembly_in_serial
virtual void complex_pre_assembly_in_serial(const model &, size_type, const model::varnamelist &, const model::varnamelist &, const model::mimlist &, model::complex_matlist &, model::complex_veclist &, model::complex_veclist &, size_type, build_version) const
Peform any pre assembly action for complex term assembly.
Definition: getfem_models.h:1551
getfem::model::interpolate_transformation_exists
bool interpolate_transformation_exists(const std::string &name) const
Tests if name corresponds to an interpolate transformation.
Definition: getfem_models.h:1129
getfem::add_normal_source_term_brick
size_type APIDECL add_normal_source_term_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataexpr, size_type region)
Add a source term on the variable varname on a boundary region.
Definition: getfem_models.cc:4279
getfem::add_basic_d_on_dt_brick
size_type APIDECL add_basic_d_on_dt_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataname_dt, const std::string &dataname_rho=std::string(), size_type region=size_type(-1))
Basic d/dt brick ( ).
Definition: getfem_models.cc:6781
getfem::model::add_initialized_fixed_size_data
void add_initialized_fixed_size_data(const std::string &name, const VECT &v)
Add a fixed size data (assumed to be a vector) to the model and initialized with v.
Definition: getfem_models.h:743
getfem_generic_assembly.h
A language for generic assembly of pde boundary value problems.
getfem::compute_isotropic_linearized_Von_Mises_pstress
void APIDECL compute_isotropic_linearized_Von_Mises_pstress(model &md, const std::string &varname, const std::string &data_E, const std::string &data_nu, const mesh_fem &mf_vm, model_real_plain_vector &VM)
Compute the Von-Mises stress of a displacement field for isotropic linearized elasticity in 3D or in ...
Definition: getfem_models.cc:6271
getfem::add_pointwise_constraints_with_given_multipliers
size_type APIDECL add_pointwise_constraints_with_given_multipliers(model &md, const std::string &varname, const std::string &multname, const std::string &dataname_pt, const std::string &dataname_unitv=std::string(), const std::string &dataname_val=std::string())
Add some pointwise constraints on the variable varname using a given multiplier multname.
Definition: getfem_models.cc:5384
getfem::model::touch_brick
void touch_brick(size_type ib)
Force the re-computation of a brick for the next assembly.
Definition: getfem_models.h:1031
getfem::velocity_update_for_Newmark_scheme
void APIDECL velocity_update_for_Newmark_scheme(model &md, size_type id2dt2b, const std::string &U, const std::string &V, const std::string &pdt, const std::string &ptwobeta, const std::string &pgamma)
Function which udpate the velocity $v^{n+1}$ after the computation of the displacement $u^{n+1}$ and ...
Definition: getfem_models.cc:7096
getfem::model::add_interpolate_transformation
void add_interpolate_transformation(const std::string &name, pinterpolate_transformation ptrans)
Add an interpolate transformation to the model to be used with the generic assembly.
Definition: getfem_models.h:1106
getfem::add_lumped_mass_for_first_order_brick
size_type APIDECL add_lumped_mass_for_first_order_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataexpr_rho=std::string(), size_type region=size_type(-1))
Lumped mass brick for first order.
Definition: getfem_models.cc:6616
getfem::context_dependencies
Deal with interdependencies of objects.
Definition: getfem_context.h:81
getfem_im_data.h
Provides indexing of integration points for mesh_im.
getfem::model::has_internal_variables
bool has_internal_variables() const
Return true if the model has at least one internal variable.
Definition: getfem_models.h:572
getfem::model::complex_tangent_matrix
const model_complex_sparse_matrix & complex_tangent_matrix() const
Gives the access to the tangent matrix.
Definition: getfem_models.h:927
getfem::model::get_active_bricks
const dal::bit_vector & get_active_bricks() const
Return the model brick ids.
Definition: getfem_models.h:1026
getfem::model::update_from_context
void update_from_context() const
this function has to be defined and should update the object when the context is modified.
Definition: getfem_models.h:506
gmm::rsvector
sparse vector built upon std::vector.
Definition: gmm_def.h:488
getfem::model::is_coercive
bool is_coercive() const
Return true if all the model terms do not affect the coercivity of the whole tangent system.
Definition: getfem_models.h:569
getfem::virtual_time_scheme
The time integration scheme object provides the necessary methods for the model object to apply a tim...
Definition: getfem_models.h:1228
getfem::model::add_secondary_domain
void add_secondary_domain(const std::string &name, psecondary_domain ptrans)
Add a secondary domain to the model to be used with the generic assembly.
Definition: getfem_models.h:1160
getfem::model::macro_dictionary
const ga_macro_dictionary & macro_dictionary() const
Dictonnary of user defined macros.
Definition: getfem_models.h:890
getfem::model::set_real_rhs
model_real_plain_vector & set_real_rhs(bool with_internal=false) const
Gives write access to the right hand side of the tangent linear system.
Definition: getfem_models.h:944
getfem::model::real_rhs
const model_real_plain_vector & real_rhs(bool with_internal=false) const
Gives access to the right hand side of the tangent linear system.
Definition: getfem_models.h:935
getfem::model::add_initialized_scalar_data
void add_initialized_scalar_data(const std::string &name, T e)
Add a scalar data (i.e.
Definition: getfem_models.h:782
getfem::model::elementary_transformation_exists
bool elementary_transformation_exists(const std::string &name) const
Tests if name corresponds to an elementary transformation.
Definition: getfem_models.h:1153
getfem::add_isotropic_linearized_elasticity_pstress_brick
size_type APIDECL add_isotropic_linearized_elasticity_pstress_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &data_E, const std::string &data_nu, size_type region)
Linear elasticity brick ( ).
Definition: getfem_models.cc:6171
getfem::virtual_brick::asm_real_tangent_terms
virtual void asm_real_tangent_terms(const model &, size_type, const model::varnamelist &, const model::varnamelist &, const model::mimlist &, model::real_matlist &, model::real_veclist &, model::real_veclist &, size_type, build_version) const
Assembly of bricks real tangent terms.
Definition: getfem_models.h:1493
getfem::add_Dirichlet_condition_with_penalization
size_type APIDECL add_Dirichlet_condition_with_penalization(model &md, const mesh_im &mim, const std::string &varname, scalar_type penalization_coeff, size_type region, const std::string &dataname=std::string(), const mesh_fem *mf_mult=0)
Add a Dirichlet condition on the variable varname and the mesh region region.
Definition: getfem_models.cc:4710
getfem::model::set_complex_rhs
model_complex_plain_vector & set_complex_rhs() const
Gives write access to the right hand side of the tangent linear system.
Definition: getfem_models.h:990
getfem::add_normal_Dirichlet_condition_with_Nitsche_method
size_type APIDECL add_normal_Dirichlet_condition_with_Nitsche_method(model &md, const mesh_im &mim, const std::string &varname, const std::string &Neumannterm, const std::string &datagamma0, size_type region, scalar_type theta=scalar_type(0), const std::string &datag=std::string())
Add a Dirichlet condition on the normal component of the variable varname and the mesh region region.
Definition: getfem_models.cc:5074
dal::static_stored_object
base class for static stored objects
Definition: dal_static_stored_objects.h:206
getfem::im_data
im_data provides indexing to the integration points of a mesh im object.
Definition: getfem_im_data.h:69
getfem::model::secondary_domain_exists
bool secondary_domain_exists(const std::string &name) const
Tests if name corresponds to an interpolate transformation.
Definition: getfem_models.h:1179
getfem::add_generalized_Dirichlet_condition_with_Nitsche_method
size_type APIDECL add_generalized_Dirichlet_condition_with_Nitsche_method(model &md, const mesh_im &mim, const std::string &varname, const std::string &Neumannterm, const std::string &datagamma0, size_type region, scalar_type theta, const std::string &datag, const std::string &dataH)
Add a Dirichlet condition on the variable varname and the mesh region region.
Definition: getfem_models.cc:5107
getfem::add_twodomain_source_term
size_type APIDECL add_twodomain_source_term(model &md, const mesh_im &mim, const std::string &expr, size_type region, const std::string &secondary_domain, const std::string &brickname=std::string(), const std::string &directvarname=std::string(), const std::string &directdataname=std::string(), bool return_if_nonlin=false)
Adds a source term given by a weak form language expression like add_source_term function but for an ...
Definition: getfem_models.cc:3434
getfem::virtual_brick::real_pre_assembly_in_serial
virtual void real_pre_assembly_in_serial(const model &, size_type, const model::varnamelist &, const model::varnamelist &, const model::mimlist &, model::real_matlist &, model::real_veclist &, model::real_veclist &, size_type, build_version) const
Peform any pre assembly action for real term assembly.
Definition: getfem_models.h:1537
getfem::virtual_brick
The virtual brick has to be derived to describe real model bricks.
Definition: getfem_models.h:1440
getfem::add_normal_Dirichlet_condition_with_penalization
size_type APIDECL add_normal_Dirichlet_condition_with_penalization(model &md, const mesh_im &mim, const std::string &varname, scalar_type penalization_coeff, size_type region, const std::string &dataname=std::string(), const mesh_fem *mf_mult=0)
Add a Dirichlet condition to the normal component of the vector (or tensor) valued variable varname a...
Definition: getfem_models.cc:4764
getfem::model::interpolate_transformation
pinterpolate_transformation interpolate_transformation(const std::string &name) const
Get a pointer to the interpolate transformation name.
Definition: getfem_models.h:1120
getfem::model::internal_solution
const model_real_plain_vector & internal_solution() const
Gives access to the partial solution for condensed internal variables.
Definition: getfem_models.h:952
getfem::add_source_term_brick
size_type APIDECL add_source_term_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataexpr, size_type region=size_type(-1), const std::string &directdataname=std::string())
Add a source term on the variable varname.
Definition: getfem_models.cc:4127
getfem::virtual_brick::asm_complex_tangent_terms
virtual void asm_complex_tangent_terms(const model &, size_type, const model::varnamelist &, const model::varnamelist &, const model::mimlist &, model::complex_matlist &, model::complex_veclist &, model::complex_veclist &, size_type, build_version) const
Assembly of bricks complex tangent terms.
Definition: getfem_models.h:1519
getfem::add_isotropic_linearized_elasticity_brick
size_type APIDECL add_isotropic_linearized_elasticity_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataname_lambda, const std::string &dataname_mu, size_type region=size_type(-1), const std::string &dataname_preconstraint=std::string())
Linear elasticity brick ( ).
Definition: getfem_models.cc:6106
getfem::add_pointwise_constraints_with_multipliers
size_type APIDECL add_pointwise_constraints_with_multipliers(model &md, const std::string &varname, const std::string &dataname_pt, const std::string &dataname_unitv=std::string(), const std::string &dataname_val=std::string())
Add some pointwise constraints on the variable varname using multiplier.
Definition: getfem_models.cc:5400
getfem::model::is_symmetric
bool is_symmetric() const
Return true if all the model terms do not affect the coercivity of the whole tangent system.
Definition: getfem_models.h:580
getfem::model::real_tangent_matrix
const model_real_sparse_matrix & real_tangent_matrix(bool internal=false) const
Gives the access to the tangent matrix.
Definition: getfem_models.h:920
getfem::model::secondary_domain
psecondary_domain secondary_domain(const std::string &name) const
Get a pointer to the interpolate transformation name.
Definition: getfem_models.h:1170
getfem::virtual_brick::real_post_assembly_in_serial
virtual void real_post_assembly_in_serial(const model &, size_type, const model::varnamelist &, const model::varnamelist &, const model::mimlist &, model::real_matlist &, model::real_veclist &, model::real_veclist &, size_type, build_version) const
Peform any post assembly action for real terms.
Definition: getfem_models.h:1565
getfem::add_basic_d2_on_dt2_brick
size_type APIDECL add_basic_d2_on_dt2_brick(model &md, const mesh_im &mim, const std::string &varnameU, const std::string &datanameV, const std::string &dataname_dt, const std::string &dataname_alpha, const std::string &dataname_rho=std::string(), size_type region=size_type(-1))
Basic d2/dt2 brick ( ).
Definition: getfem_models.cc:6956
getfem::add_Laplacian_brick
size_type APIDECL add_Laplacian_brick(model &md, const mesh_im &mim, const std::string &varname, size_type region=size_type(-1))
Add a Laplacian term on the variable varname (in fact with a minus : :math:-\text{div}(\nabla u)).
Definition: getfem_models.cc:3921
getfem::pbrick
std::shared_ptr< const virtual_brick > pbrick
type of pointer on a brick
Definition: getfem_models.h:49
getfem::add_Dirichlet_condition_with_simplification
size_type APIDECL add_Dirichlet_condition_with_simplification(model &md, const std::string &varname, size_type region, const std::string &dataname=std::string())
Add a (simple) Dirichlet condition on the variable varname and the mesh region region.
Definition: getfem_models.cc:5023
getfem::add_Fourier_Robin_brick
size_type APIDECL add_Fourier_Robin_brick(model &md, const mesh_im &mim, const std::string &varname, const std::string &dataexpr, size_type region)
Add a Fourier-Robin brick to the model.
Definition: getfem_models.cc:5637
getfem::is_old
bool is_old(const std::string &name)
Does the variable have Old_ prefix.
Definition: getfem_models.cc:218
getfem::mesh_fem::nb_dof
virtual size_type nb_dof() const
Return the total number of degrees of freedom.
Definition: getfem_mesh_fem.h:562
getfem_assembling.h
Miscelleanous assembly routines for common terms. Use the low-level generic assembly....
getfem::add_source_term
size_type APIDECL add_source_term(model &md, const mesh_im &mim, const std::string &expr, size_type region=size_type(-1), const std::string &brickname=std::string(), const std::string &directvarname=std::string(), const std::string &directdataname=std::string(), bool return_if_nonlin=false)
Add a source term given by the assembly string expr which will be assembled in region region and with...
Definition: getfem_models.cc:3426
getfem::no_old_prefix_name
std::string no_old_prefix_name(const std::string &name)
Strip the variable name from prefix Old_ if it has one.
Definition: getfem_models.cc:222