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