GetFEM  5.4.3
getfem_generic_assembly_tree.cc
1 /*===========================================================================
2 
3  Copyright (C) 2013-2020 Yves Renard
4 
5  This file is a part of GetFEM
6 
7  GetFEM is free software; you can redistribute it and/or modify it
8  under the terms of the GNU Lesser General Public License as published
9  by the Free Software Foundation; either version 3 of the License, or
10  (at your option) any later version along with the GCC Runtime Library
11  Exception either version 3.1 or (at your option) any later version.
12  This program is distributed in the hope that it will be useful, but
13  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15  License and GCC Runtime Library Exception for more details.
16  You should have received a copy of the GNU Lesser General Public License
17  along with this program; if not, write to the Free Software Foundation,
18  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
19 
20 ===========================================================================*/
21 
22 #include "getfem/getfem_generic_assembly_tree.h"
23 #include "getfem/getfem_generic_assembly_functions_and_operators.h"
24 #include "getfem/getfem_generic_assembly_compile_and_exec.h"
25 
26 
27 namespace getfem {
28 
29  //=========================================================================
30  // Lexical analysis for the generic assembly language
31  //=========================================================================
32 
33  static GA_TOKEN_TYPE ga_char_type[256];
34  static int ga_operator_priorities[GA_NB_TOKEN_TYPE];
35 
36  // Initialize ga_char_type and ga_operator_priorities arrays
37  static bool init_ga_char_type() {
38  for (int i = 0; i < 256; ++i) ga_char_type[i] = GA_INVALID;
39  ga_char_type['+'] = GA_PLUS; ga_char_type['-'] = GA_MINUS;
40  ga_char_type['*'] = GA_MULT; ga_char_type['/'] = GA_DIV;
41  ga_char_type[':'] = GA_COLON; ga_char_type['\''] = GA_QUOTE;
42  ga_char_type['.'] = GA_DOT; ga_char_type['@'] = GA_TMULT;
43  ga_char_type[','] = GA_COMMA; ga_char_type[';'] = GA_SEMICOLON;
44  ga_char_type['('] = GA_LPAR; ga_char_type[')'] = GA_RPAR;
45  ga_char_type['['] = GA_LBRACKET; ga_char_type[']'] = GA_RBRACKET;
46  ga_char_type['_'] = GA_NAME; ga_char_type['='] = GA_COLON_EQ;
47  for (unsigned i = 'a'; i <= 'z'; ++i) ga_char_type[i] = GA_NAME;
48  for (unsigned i = 'A'; i <= 'Z'; ++i) ga_char_type[i] = GA_NAME;
49  for (unsigned i = '0'; i <= '9'; ++i) ga_char_type[i] = GA_SCALAR;
50 
51  for (unsigned i=0; i < GA_NB_TOKEN_TYPE; ++i) ga_operator_priorities[i] = 0;
52  ga_operator_priorities[GA_PLUS] = 1;
53  ga_operator_priorities[GA_MINUS] = 1;
54  ga_operator_priorities[GA_MULT] = 2;
55  ga_operator_priorities[GA_DIV] = 2;
56  ga_operator_priorities[GA_COLON] = 2;
57  ga_operator_priorities[GA_DOT] = 2;
58  ga_operator_priorities[GA_DOTMULT] = 2;
59  ga_operator_priorities[GA_DOTDIV] = 2;
60  ga_operator_priorities[GA_TMULT] = 2;
61  ga_operator_priorities[GA_QUOTE] = 3;
62  ga_operator_priorities[GA_UNARY_MINUS] = 3;
63  ga_operator_priorities[GA_SYM] = 4;
64  ga_operator_priorities[GA_SKEW] = 4;
65  ga_operator_priorities[GA_TRACE] = 4;
66  ga_operator_priorities[GA_DEVIATOR] = 4;
67  ga_operator_priorities[GA_PRINT] = 4;
68 
69  return true;
70  }
71 
72  static bool ga_initialized = init_ga_char_type();
73 
74  // Get the next token in the string at position 'pos' end return its type
75  static GA_TOKEN_TYPE ga_get_token(const std::string &expr,
76  size_type &pos,
77  size_type &token_pos,
78  size_type &token_length) {
79  bool fdot = false, fE = false;
80  GMM_ASSERT1(ga_initialized, "Internal error");
81 
82  // Ignore white spaces
83  while (expr[pos] == ' ' && pos < expr.size()) ++pos;
84  token_pos = pos;
85  token_length = 0;
86 
87  // Detecting end of expression
88  if (pos >= expr.size()) return GA_END;
89 
90  // Treating the different cases (Operation, name or number)
91  GA_TOKEN_TYPE type = ga_char_type[unsigned(expr[pos])];
92  ++pos; ++token_length;
93  if (type == GA_DOT) {
94  if (pos >= expr.size()) return type;
95  if (expr[pos] == '*') { ++pos; ++token_length; return GA_DOTMULT; }
96  if (expr[pos] == '/') { ++pos; ++token_length; return GA_DOTDIV; }
97  if (ga_char_type[unsigned(expr[pos])] != GA_SCALAR) return type;
98  fdot = true; type = GA_SCALAR;
99  }
100  switch (type) {
101  case GA_SCALAR:
102  while (pos < expr.size()) {
103  GA_TOKEN_TYPE ctype = ga_char_type[unsigned(expr[pos])];
104  switch (ctype) {
105  case GA_DOT:
106  if (fdot) return type;
107  fdot = true; ++pos; ++token_length;
108  break;
109  case GA_NAME:
110  if (fE || (expr[pos] != 'E' && expr[pos] != 'e')) return type;
111  fE = true; fdot = true; ++pos; ++token_length;
112  if (pos < expr.size()) {
113  if (expr[pos] == '+' || expr[pos] == '-')
114  { ++pos; ++token_length; }
115  }
116  if (pos >= expr.size()
117  || ga_char_type[unsigned(expr[pos])] != GA_SCALAR)
118  return GA_INVALID;
119  break;
120  case GA_SCALAR:
121  ++pos; ++token_length; break;
122  default:
123  return type;
124  }
125  }
126  return type;
127  case GA_NAME:
128  while (pos < expr.size()) {
129  GA_TOKEN_TYPE ctype = ga_char_type[unsigned(expr[pos])];
130  if (ctype != GA_SCALAR && ctype != GA_NAME) break;
131  ++pos; ++token_length;
132  }
133  if (expr.compare(token_pos, token_length, "Sym") == 0)
134  return GA_SYM;
135  if (expr.compare(token_pos, token_length, "Def") == 0)
136  return GA_DEF;
137  if (expr.compare(token_pos, token_length, "Skew") == 0)
138  return GA_SKEW;
139  if (expr.compare(token_pos, token_length, "Trace") == 0)
140  return GA_TRACE;
141  if (expr.compare(token_pos, token_length, "Deviator") == 0)
142  return GA_DEVIATOR;
143  if (expr.compare(token_pos, token_length, "Interpolate") == 0)
144  return GA_INTERPOLATE;
145  if (expr.compare(token_pos, token_length, "Interpolate_derivative") == 0)
146  return GA_INTERPOLATE_DERIVATIVE;
147  if (expr.compare(token_pos, token_length, "Interpolate_filter") == 0)
148  return GA_INTERPOLATE_FILTER;
149  if (expr.compare(token_pos, token_length,
150  "Elementary_transformation") == 0)
151  return GA_ELEMENTARY;
152  if (expr.compare(token_pos, token_length, "Secondary_domain") == 0 ||
153  expr.compare(token_pos, token_length, "Secondary_Domain") == 0)
154  return GA_SECONDARY_DOMAIN;
155  if (expr.compare(token_pos, token_length, "Xfem_plus") == 0)
156  return GA_XFEM_PLUS;
157  if (expr.compare(token_pos, token_length, "Xfem_minus") == 0)
158  return GA_XFEM_MINUS;
159  if (expr.compare(token_pos, token_length, "Print") == 0)
160  return GA_PRINT;
161  return type;
162  case GA_COMMA:
163  if (pos < expr.size() &&
164  ga_char_type[unsigned(expr[pos])] == GA_COMMA) {
165  ++pos; return GA_DCOMMA;
166  }
167  return type;
168  case GA_SEMICOLON:
169  if (pos < expr.size() &&
170  ga_char_type[unsigned(expr[pos])] == GA_SEMICOLON) {
171  ++pos; return GA_DSEMICOLON;
172  }
173  return type;
174  case GA_COLON:
175  if (pos < expr.size() &&
176  ga_char_type[unsigned(expr[pos])] == GA_COLON_EQ) {
177  ++pos; return GA_COLON_EQ;
178  }
179  return type;
180  case GA_COLON_EQ:
181  return GA_INVALID;
182  default: return type;
183  }
184  }
185 
186  //=========================================================================
187  // Error handling
188  //=========================================================================
189 
190  void ga_throw_error_msg(pstring expr, size_type pos,
191  const std::string &msg) {
192  int length_before = 70, length_after = 70;
193  if (expr && expr->size()) {
194  int first = std::max(0, int(pos)-length_before);
195  int last = std::min(int(pos)+length_after, int(expr->size()));
196  if (last - first < length_before+length_after)
197  first = std::max(0, int(pos)-length_before
198  -(length_before+length_after-last+first));
199  if (last - first < length_before+length_after)
200  last = std::min(int(pos)+length_after
201  +(length_before+length_after-last+first),
202  int(expr->size()));
203  if (first > 0) cerr << "...";
204  cerr << expr->substr(first, last-first);
205  if (last < int(expr->size())) cerr << "...";
206  cerr << endl;
207  if (first > 0) cerr << " ";
208  if (int(pos) > first)
209  cerr << std::setfill ('-') << std::setw(int(pos)-first) << '-'
210  << std::setfill (' ');
211  cerr << "^" << endl;
212  }
213  cerr << msg << endl;
214  }
215 
216  //=========================================================================
217  // Tree structure
218  //=========================================================================
219 
220  void ga_tree_node::mult_test(const pga_tree_node n0, const pga_tree_node n1) {
221 
222  size_type test0 = n0->test_function_type, test1 = n1->test_function_type;
223  if (test0 && test1 && (test0 == test1 || test0 >= 3 || test1 >= 3))
224  ga_throw_error(expr, pos,
225  "Incompatibility of test functions in product.");
226  GMM_ASSERT1(test0 != size_type(-1) && test1 != size_type(-1),
227  "internal error");
228 
229  test_function_type = test0 + test1;
230 
231  size_type st = nb_test_functions();
232  bgeot::multi_index mi(st);
233 
234  switch (test0) {
235  case 1: mi[0] = n0->t.sizes()[0]; break;
236  case 2: mi[st-1] = n0->t.sizes()[0]; break;
237  case 3: mi[0] = n0->t.sizes()[0]; mi[1] = n0->t.sizes()[1]; break;
238  }
239  switch (test1) {
240  case 1: mi[0] = n1->t.sizes()[0]; break;
241  case 2: mi[st-1] = n1->t.sizes()[0]; break;
242  case 3: mi[0] = n1->t.sizes()[0]; mi[1] = n1->t.sizes()[1]; break;
243  }
244 
245  if (n0->name_test1.size()) {
246  name_test1 = n0->name_test1; qdim1 = n0->qdim1;
247  interpolate_name_test1 = n0->interpolate_name_test1;
248  } else {
249  name_test1 = n1->name_test1; qdim1 = n1->qdim1;
250  interpolate_name_test1 = n1->interpolate_name_test1;
251  }
252 
253  if (n0->name_test2.size()) {
254  name_test2 = n0->name_test2; qdim2 = n0->qdim2;
255  interpolate_name_test2 = n0->interpolate_name_test2;
256  } else {
257  name_test2 = n1->name_test2; qdim2 = n1->qdim2;
258  interpolate_name_test2 = n1->interpolate_name_test2;
259  }
260  t.adjust_sizes(mi);
261  }
262 
263  void ga_tree::add_scalar(scalar_type val, size_type pos, pstring expr) {
264  while (current_node && current_node->node_type != GA_NODE_OP)
265  current_node = current_node->parent;
266  if (current_node) {
267  current_node->adopt_child(new ga_tree_node(val, pos, expr));
268  current_node = current_node->children.back();
269  }
270  else {
271  GMM_ASSERT1(root == nullptr, "Invalid tree operation");
272  current_node = root = new ga_tree_node(val, pos, expr);
273  root->parent = nullptr;
274  }
275  }
276 
277  void ga_tree::add_allindices(size_type pos, pstring expr) {
278  while (current_node && current_node->node_type != GA_NODE_OP)
279  current_node = current_node->parent;
280  if (current_node) {
281  current_node->adopt_child(new ga_tree_node(GA_NODE_ALLINDICES, pos,expr));
282  current_node = current_node->children.back();
283  }
284  else {
285  GMM_ASSERT1(root == nullptr, "Invalid tree operation");
286  current_node = root = new ga_tree_node(GA_NODE_ALLINDICES, pos, expr);
287  root->parent = nullptr;
288  }
289  }
290 
291  void ga_tree::add_name(const char *name, size_type length, size_type pos,
292  pstring expr) {
293  while (current_node && current_node->node_type != GA_NODE_OP)
294  current_node = current_node->parent;
295  if (current_node) {
296  current_node->adopt_child(new ga_tree_node(name, length, pos, expr));
297  current_node = current_node->children.back();
298  }
299  else {
300  GMM_ASSERT1(root == nullptr, "Invalid tree operation");
301  current_node = root = new ga_tree_node(name, length, pos, expr);
302  root->parent = nullptr;
303  }
304  }
305 
306  void ga_tree::add_sub_tree(ga_tree &sub_tree) {
307  if (current_node &&
308  (current_node->node_type == GA_NODE_PARAMS ||
309  current_node->node_type == GA_NODE_INTERPOLATE_FILTER ||
310  current_node->node_type == GA_NODE_C_MATRIX)) {
311  GMM_ASSERT1(sub_tree.root, "Invalid tree operation");
312  current_node->adopt_child(sub_tree.root);
313  } else {
314  GMM_ASSERT1(sub_tree.root, "Invalid tree operation");
315  while (current_node && current_node->node_type != GA_NODE_OP)
316  current_node = current_node->parent;
317  if (current_node) {
318  current_node->adopt_child(sub_tree.root);
319  current_node = sub_tree.root;
320  }
321  else {
322  GMM_ASSERT1(root == nullptr, "Invalid tree operation");
323  current_node = root = sub_tree.root;
324  root->parent = nullptr;
325  }
326  }
327  sub_tree.root = sub_tree.current_node = nullptr;
328  }
329 
330  void ga_tree::add_params(size_type pos, pstring expr) {
331  GMM_ASSERT1(current_node, "internal error");
332  while (current_node && current_node->parent &&
333  current_node->parent->node_type == GA_NODE_OP &&
334  ga_operator_priorities[current_node->parent->op_type] >= 4)
335  current_node = current_node->parent;
336  pga_tree_node new_node = new ga_tree_node(GA_NODE_PARAMS, pos, expr);
337  new_node->parent = current_node->parent;
338  if (current_node->parent)
339  current_node->parent->replace_child(current_node, new_node);
340  else
341  root = new_node;
342  new_node->adopt_child(current_node);
343  current_node = new_node;
344  }
345 
346  void ga_tree::add_matrix(size_type pos, pstring expr) {
347  while (current_node && current_node->node_type != GA_NODE_OP)
348  current_node = current_node->parent;
349  if (current_node) {
350  current_node->adopt_child(new ga_tree_node(GA_NODE_C_MATRIX, pos, expr));
351  current_node = current_node->children.back();
352  }
353  else {
354  GMM_ASSERT1(root == nullptr, "Invalid tree operation");
355  current_node = root = new ga_tree_node(GA_NODE_C_MATRIX, pos, expr);
356  root->parent = nullptr;
357  }
358  current_node->nbc1 = current_node->nbc2 = current_node->nbc3 = 0;
359  }
360 
361  void ga_tree::add_op(GA_TOKEN_TYPE op_type, size_type pos,
362  pstring expr) {
363  while (current_node && current_node->parent &&
364  current_node->parent->node_type == GA_NODE_OP &&
365  ga_operator_priorities[current_node->parent->op_type]
366  >= ga_operator_priorities[op_type])
367  current_node = current_node->parent;
368  pga_tree_node new_node = new ga_tree_node(op_type, pos, expr);
369  if (current_node) {
370  if (op_type == GA_UNARY_MINUS
371  || op_type == GA_SYM || op_type == GA_SKEW
372  || op_type == GA_TRACE || op_type == GA_DEVIATOR
373  || op_type == GA_PRINT) {
374  current_node->adopt_child(new_node);
375  } else {
376  new_node->parent = current_node->parent;
377  if (current_node->parent)
378  current_node->parent->replace_child(current_node, new_node);
379  else
380  root = new_node;
381  new_node->adopt_child(current_node);
382  }
383  } else {
384  if (root) new_node->adopt_child(root);
385  root = new_node;
386  root->parent = nullptr;
387  }
388  current_node = new_node;
389  }
390 
391  void ga_tree::clear_node_rec(pga_tree_node pnode) {
392  if (pnode) {
393  for (pga_tree_node &child : pnode->children)
394  clear_node_rec(child);
395  delete pnode;
396  current_node = nullptr;
397  }
398  }
399 
400  void ga_tree::clear_node(pga_tree_node pnode) {
401  if (pnode) {
402  pga_tree_node parent = pnode->parent;
403  if (parent) { // keep all siblings of pnode
404  size_type j = 0;
405  for (pga_tree_node &sibling : parent->children)
406  if (sibling != pnode)
407  parent->children[j++] = sibling;
408  parent->children.resize(j);
409  } else
410  root = nullptr;
411  }
412  clear_node_rec(pnode);
413  }
414 
415  void ga_tree::clear_children(pga_tree_node pnode) {
416  for (pga_tree_node &child : pnode->children)
417  clear_node_rec(child);
418  pnode->children.resize(0);
419  }
420 
421  void ga_tree::replace_node_by_child(pga_tree_node pnode, size_type i) {
422  GMM_ASSERT1(i < pnode->children.size(), "Internal error");
423  pga_tree_node child = pnode->children[i];
424  child->parent = pnode->parent;
425  if (pnode->parent)
426  pnode->parent->replace_child(pnode, child);
427  else
428  root = child;
429  current_node = nullptr;
430  for (pga_tree_node &sibling : pnode->children)
431  if (sibling != child) clear_node_rec(sibling);
432  delete pnode;
433  }
434 
435  // This function allocates a new node "pnode_new", copies the content of "pnode"
436  // into the newly allocated node (including deep copies of all of its children)
437  void ga_tree::copy_node(pga_tree_node pnode,
438  pga_tree_node &pnode_new) {
439  GMM_ASSERT1(pnode_new == nullptr, "Internal error");
440  pnode_new = new ga_tree_node();
441  *pnode_new = *pnode;
442  pnode_new->parent = nullptr;
443  for (size_type j = 0; j < pnode_new->children.size(); ++j) {
444  pnode_new->children[j] = nullptr;
445  copy_node(pnode->children[j], pnode_new->children[j]);
446  pnode_new->accept_child(j);
447  }
448  }
449 
450  void ga_tree::duplicate_with_operation(pga_tree_node pnode,
451  GA_TOKEN_TYPE op_type) {
452  pga_tree_node newop = new ga_tree_node(op_type, pnode->pos, pnode->expr);
453  newop->children.resize(2, nullptr);
454  newop->children[0] = pnode;
455  newop->pos = pnode->pos; newop->expr = pnode->expr;
456  newop->parent = pnode->parent;
457  if (pnode->parent)
458  pnode->parent->replace_child(pnode, newop);
459  else
460  root = newop;
461  pnode->parent = newop;
462  copy_node(pnode, newop->children[1]);
463  newop->accept_child(1);
464  }
465 
466  void ga_tree::add_child(pga_tree_node pnode, GA_NODE_TYPE node_type) {
467  pga_tree_node newnode=new ga_tree_node();
468  newnode->pos = pnode->pos; newnode->expr = pnode->expr;
469  newnode->node_type = node_type; pnode->adopt_child(newnode);
470  }
471 
472  void ga_tree::insert_node(pga_tree_node pnode, GA_NODE_TYPE node_type) {
473  pga_tree_node newnode = new ga_tree_node();
474  newnode->node_type = node_type;
475  newnode->parent = pnode->parent;
476  newnode->pos = pnode->pos; newnode->expr = pnode->expr;
477  if (pnode->parent)
478  pnode->parent->replace_child(pnode, newnode);
479  else
480  root = newnode;
481  newnode->adopt_child(pnode);
482  }
483 
484  bool sub_tree_are_equal
485  (const pga_tree_node pnode1, const pga_tree_node pnode2,
486  const ga_workspace &workspace, int version) {
487 
488  size_type ntype1 = pnode1->node_type;
489  if (ntype1 == GA_NODE_ZERO) ntype1 = GA_NODE_CONSTANT;
490  size_type ntype2 = pnode2->node_type;
491  if (ntype2 == GA_NODE_ZERO) ntype2 = GA_NODE_CONSTANT;
492  if (ntype1 != ntype2) return false;
493  if (pnode1->children.size() != pnode2->children.size()) return false;
494 
495  switch(ntype1) {
496  case GA_NODE_OP:
497  if (pnode1->op_type != pnode2->op_type) return false;
498  if (pnode1->symmetric_op != pnode2->symmetric_op) return false;
499  break;
500  case GA_NODE_OPERATOR:
501  if (pnode1->der1 != pnode2->der1 || pnode1->der2 != pnode2->der2)
502  return false;
503  if (pnode1->name.compare(pnode2->name)) return false;
504  break;
505  case GA_NODE_PREDEF_FUNC: case GA_NODE_SPEC_FUNC:
506  if (pnode1->name.compare(pnode2->name)) return false;
507  break;
508  case GA_NODE_CONSTANT: case GA_NODE_ZERO:
509  if (pnode1->tensor().size() != pnode2->tensor().size()) return false;
510 
511  switch(version) {
512  case 0: case 1:
513  if (pnode1->test_function_type != pnode2->test_function_type)
514  return false;
515  if ((pnode1->test_function_type & 1) &&
516  pnode1->name_test1.compare(pnode2->name_test1) != 0)
517  return false;
518  if ((pnode1->test_function_type & 2) &&
519  pnode1->name_test2.compare(pnode2->name_test2) != 0)
520  return false;
521  break;
522  case 2:
523  if ((pnode1->test_function_type == 1 &&
524  pnode2->test_function_type == 1) ||
525  (pnode1->test_function_type == 2 &&
526  pnode2->test_function_type == 2))
527  return false;
528  if ((pnode1->test_function_type & 1) &&
529  pnode1->name_test1.compare(pnode2->name_test2) != 0)
530  return false;
531  if ((pnode1->test_function_type & 2) &&
532  pnode1->name_test2.compare(pnode2->name_test1) != 0)
533  return false;
534  break;
535  }
536  if (pnode1->tensor().size() != 1 &&
537  pnode1->t.sizes().size() != pnode2->t.sizes().size()) return false;
538  for (size_type i = 0; i < pnode1->t.sizes().size(); ++i)
539  if (pnode1->t.sizes()[i] != pnode2->t.sizes()[i]) return false;
540  for (size_type i = 0; i < pnode1->tensor().size(); ++i)
541  if (gmm::abs(pnode1->tensor()[i] - pnode2->tensor()[i]) > 1E-25)
542  return false;
543  break;
544  case GA_NODE_C_MATRIX:
545  if (pnode1->children.size() != pnode2->children.size() ||
546  pnode1->nb_test_functions() != pnode2->nb_test_functions() ||
547  pnode1->t.sizes().size() != pnode2->t.sizes().size())
548  return false;
549  for (size_type i=pnode1->nb_test_functions();
550  i < pnode1->t.sizes().size(); ++i)
551  if (pnode1->t.sizes()[i] != pnode2->t.sizes()[i]) return false;
552  if (pnode1->nbc1 != pnode2->nbc1) return false;
553  break;
554  case GA_NODE_INTERPOLATE_FILTER:
555  if (pnode1->interpolate_name.compare(pnode2->interpolate_name) ||
556  pnode1->nbc1 != pnode2->nbc1)
557  return false;
558  break;
559  case GA_NODE_INTERPOLATE_X:
560  case GA_NODE_SECONDARY_DOMAIN_X:
561  case GA_NODE_INTERPOLATE_ELT_K: case GA_NODE_INTERPOLATE_ELT_B:
562  case GA_NODE_INTERPOLATE_NORMAL:
563  case GA_NODE_SECONDARY_DOMAIN_NORMAL:
564  if (pnode1->interpolate_name.compare(pnode2->interpolate_name))
565  return false;
566  break;
567  case GA_NODE_INTERPOLATE_DERIVATIVE:
568  if (pnode1->interpolate_name_der.compare(pnode2->interpolate_name_der))
569  return false;
570  if (pnode1->interpolate_name.compare(pnode2->interpolate_name) ||
571  pnode1->elementary_name.compare(pnode2->elementary_name))
572  return false;
573  {
574  const mesh_fem *mf1 = workspace.associated_mf(pnode1->name);
575  const mesh_fem *mf2 = workspace.associated_mf(pnode2->name);
576  switch (version) {
577  case 0:
578  if (pnode1->name.compare(pnode2->name) ||
579  pnode1->test_function_type != pnode2->test_function_type)
580  return false;
581  break;
582  case 1:
583  if (mf1 != mf2 ||
584  workspace.qdim(pnode1->name) != workspace.qdim(pnode2->name) ||
585  pnode1->test_function_type != pnode2->test_function_type)
586  return false;
587  break;
588  case 2:
589  if (mf1 != mf2 ||
590  workspace.qdim(pnode1->name) != workspace.qdim(pnode2->name) ||
591  pnode1->test_function_type == pnode2->test_function_type)
592  return false;
593  break;
594  }
595  }
596  break;
597  case GA_NODE_INTERPOLATE_VAL_TEST: case GA_NODE_INTERPOLATE_GRAD_TEST:
598  case GA_NODE_INTERPOLATE_HESS_TEST: case GA_NODE_INTERPOLATE_DIVERG_TEST:
599  case GA_NODE_ELEMENTARY_VAL_TEST: case GA_NODE_ELEMENTARY_GRAD_TEST:
600  case GA_NODE_ELEMENTARY_HESS_TEST: case GA_NODE_ELEMENTARY_DIVERG_TEST:
601  case GA_NODE_SECONDARY_DOMAIN_VAL_TEST:
602  case GA_NODE_SECONDARY_DOMAIN_GRAD_TEST:
603  case GA_NODE_SECONDARY_DOMAIN_HESS_TEST:
604  case GA_NODE_SECONDARY_DOMAIN_DIVERG_TEST:
605  case GA_NODE_XFEM_PLUS_VAL_TEST: case GA_NODE_XFEM_PLUS_GRAD_TEST:
606  case GA_NODE_XFEM_PLUS_HESS_TEST: case GA_NODE_XFEM_PLUS_DIVERG_TEST:
607  case GA_NODE_XFEM_MINUS_VAL_TEST: case GA_NODE_XFEM_MINUS_GRAD_TEST:
608  case GA_NODE_XFEM_MINUS_HESS_TEST: case GA_NODE_XFEM_MINUS_DIVERG_TEST:
609  if (pnode1->interpolate_name.compare(pnode2->interpolate_name) ||
610  pnode1->elementary_name.compare(pnode2->elementary_name))
611  return false;
612  {
613  const mesh_fem *mf1 = workspace.associated_mf(pnode1->name);
614  const mesh_fem *mf2 = workspace.associated_mf(pnode2->name);
615  switch (version) {
616  case 0:
617  if (pnode1->name.compare(pnode2->name) ||
618  pnode1->test_function_type != pnode2->test_function_type)
619  return false;
620  break;
621  case 1:
622  if (mf1 != mf2 ||
623  workspace.qdim(pnode1->name) != workspace.qdim(pnode2->name) ||
624  pnode1->test_function_type != pnode2->test_function_type)
625  return false;
626  break;
627  case 2:
628  if (mf1 != mf2 ||
629  workspace.qdim(pnode1->name) != workspace.qdim(pnode2->name) ||
630  pnode1->test_function_type == pnode2->test_function_type)
631  return false;
632  break;
633  }
634  }
635  break;
636  case GA_NODE_VAL_TEST: case GA_NODE_GRAD_TEST:
637  case GA_NODE_HESS_TEST: case GA_NODE_DIVERG_TEST:
638  {
639  const mesh_fem *mf1 = workspace.associated_mf(pnode1->name);
640  const mesh_fem *mf2 = workspace.associated_mf(pnode2->name);
641  switch (version) {
642  case 0:
643  if (pnode1->name.compare(pnode2->name) ||
644  pnode1->test_function_type != pnode2->test_function_type)
645  return false;
646  break;
647  case 1:
648  if (mf1 != mf2 ||
649  workspace.qdim(pnode1->name) != workspace.qdim(pnode2->name) ||
650  pnode1->test_function_type != pnode2->test_function_type)
651  return false;
652  break;
653  case 2:
654  if (mf1 != mf2 ||
655  workspace.qdim(pnode1->name) != workspace.qdim(pnode2->name) ||
656  pnode1->test_function_type == pnode2->test_function_type)
657  return false;
658  break;
659  }
660  }
661  break;
662  case GA_NODE_VAL: case GA_NODE_GRAD:
663  case GA_NODE_HESS: case GA_NODE_DIVERG:
664  if (pnode1->name.compare(pnode2->name)) return false;
665  break;
666  case GA_NODE_INTERPOLATE_VAL: case GA_NODE_INTERPOLATE_GRAD:
667  case GA_NODE_INTERPOLATE_HESS: case GA_NODE_INTERPOLATE_DIVERG:
668  case GA_NODE_ELEMENTARY_VAL: case GA_NODE_ELEMENTARY_GRAD:
669  case GA_NODE_ELEMENTARY_HESS: case GA_NODE_ELEMENTARY_DIVERG:
670  case GA_NODE_SECONDARY_DOMAIN_VAL: case GA_NODE_SECONDARY_DOMAIN_GRAD:
671  case GA_NODE_SECONDARY_DOMAIN_HESS: case GA_NODE_SECONDARY_DOMAIN_DIVERG:
672  case GA_NODE_XFEM_PLUS_VAL: case GA_NODE_XFEM_PLUS_GRAD:
673  case GA_NODE_XFEM_PLUS_HESS: case GA_NODE_XFEM_PLUS_DIVERG:
674  case GA_NODE_XFEM_MINUS_VAL: case GA_NODE_XFEM_MINUS_GRAD:
675  case GA_NODE_XFEM_MINUS_HESS: case GA_NODE_XFEM_MINUS_DIVERG:
676  if (pnode1->interpolate_name.compare(pnode2->interpolate_name) ||
677  pnode1->elementary_name.compare(pnode2->elementary_name) ||
678  pnode1->name.compare(pnode2->name))
679  return false;
680  break;
681  case GA_NODE_X:
682  if (pnode1->nbc1 != pnode2->nbc1) return false;
683  break;
684 
685  default:break;
686  }
687 
688  if (version && ntype1 == GA_NODE_OP && pnode1->symmetric_op) {
689  if (sub_tree_are_equal(pnode1->children[0], pnode2->children[0],
690  workspace, version) &&
691  sub_tree_are_equal(pnode1->children[1], pnode2->children[1],
692  workspace, version))
693  return true;
694  if (sub_tree_are_equal(pnode1->children[1], pnode2->children[0],
695  workspace, version) &&
696  sub_tree_are_equal(pnode1->children[0], pnode2->children[1],
697  workspace, version) )
698  return true;
699  return false;
700  } else {
701  for (size_type i = 0; i < pnode1->children.size(); ++i)
702  if (!(sub_tree_are_equal(pnode1->children[i], pnode2->children[i],
703  workspace, version)))
704  return false;
705  }
706  return true;
707  }
708 
709  static void verify_tree(const pga_tree_node pnode,
710  const pga_tree_node parent) {
711  GMM_ASSERT1(pnode->parent == parent,
712  "Invalid tree node " << pnode->node_type);
713  for (pga_tree_node &child : pnode->children)
714  verify_tree(child, pnode);
715  }
716 
717 
718  static void ga_print_constant_tensor(const pga_tree_node pnode,
719  std::ostream &str) {
720  size_type nt = pnode->nb_test_functions(); // for printing zero tensors
721  switch (pnode->tensor_order()) {
722  case 0:
723  str << (nt ? scalar_type(0) : pnode->tensor()[0]);
724  break;
725 
726  case 1:
727  str << "[";
728  for (size_type i = 0; i < pnode->tensor_proper_size(0); ++i) {
729  if (i != 0) str << ",";
730  str << (nt ? scalar_type(0) : pnode->tensor()[i]);
731  }
732  str << "]";
733  break;
734 
735  case 2: case 3: case 4:
736  {
737  size_type ii(0);
738  size_type n0 = pnode->tensor_proper_size(0);
739  size_type n1 = pnode->tensor_proper_size(1);
740  size_type n2 = ((pnode->tensor_order() > 2) ?
741  pnode->tensor_proper_size(2) : 1);
742  size_type n3 = ((pnode->tensor_order() > 3) ?
743  pnode->tensor_proper_size(3) : 1);
744  if (n3 > 1) str << "[";
745  for (size_type l = 0; l < n3; ++l) {
746  if (l != 0) str << ",";
747  if (n2 > 1) str << "[";
748  for (size_type k = 0; k < n2; ++k) {
749  if (k != 0) str << ",";
750  if (n1 > 1) str << "[";
751  for (size_type j = 0; j < n1; ++j) {
752  if (j != 0) str << ",";
753  if (n0 > 1) str << "[";
754  for (size_type i = 0; i < n0; ++i) {
755  if (i != 0) str << ",";
756  str << (nt ? scalar_type(0) : pnode->tensor()[ii++]);
757  }
758  if (n0 > 1) str << "]";
759  }
760  if (n1 > 1) str << "]";
761  }
762  if (n2 > 1) str << "]";
763  }
764  if (n3 > 1) str << "]";
765  }
766  break;
767 
768  case 5: case 6: case 7: case 8:
769  str << "Reshape([";
770  for (size_type i = 0; i < pnode->tensor_proper_size(); ++i) {
771  if (i != 0) str << ";";
772  str << (nt ? scalar_type(0) : pnode->tensor()[i]);
773  }
774  str << "]";
775  for (size_type i = 0; i < pnode->tensor_order(); ++i) {
776  if (i != 0) str << ", ";
777  str << pnode->tensor_proper_size(i);
778  }
779  str << ")";
780  break;
781 
782  default: GMM_ASSERT1(false, "Invalid tensor dimension");
783  }
784  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
785  }
786 
787  void ga_print_node(const pga_tree_node pnode, std::ostream &str) {
788  if (!pnode) return;
789  long prec = str.precision(16);
790 
791  bool is_interpolate(false), is_elementary(false), is_secondary(false);
792  bool is_xfem_plus(false), is_xfem_minus(false);
793  switch(pnode->node_type) {
794  case GA_NODE_INTERPOLATE:
795  case GA_NODE_INTERPOLATE_X:
796  case GA_NODE_INTERPOLATE_ELT_K: case GA_NODE_INTERPOLATE_ELT_B:
797  case GA_NODE_INTERPOLATE_NORMAL:
798  case GA_NODE_INTERPOLATE_VAL:
799  case GA_NODE_INTERPOLATE_GRAD:
800  case GA_NODE_INTERPOLATE_HESS:
801  case GA_NODE_INTERPOLATE_DIVERG:
802  case GA_NODE_INTERPOLATE_VAL_TEST:
803  case GA_NODE_INTERPOLATE_GRAD_TEST:
804  case GA_NODE_INTERPOLATE_HESS_TEST:
805  case GA_NODE_INTERPOLATE_DIVERG_TEST:
806  str << "Interpolate(";
807  is_interpolate = true;
808  break;
809  case GA_NODE_ELEMENTARY:
810  case GA_NODE_ELEMENTARY_VAL:
811  case GA_NODE_ELEMENTARY_GRAD:
812  case GA_NODE_ELEMENTARY_HESS:
813  case GA_NODE_ELEMENTARY_DIVERG:
814  case GA_NODE_ELEMENTARY_VAL_TEST:
815  case GA_NODE_ELEMENTARY_GRAD_TEST:
816  case GA_NODE_ELEMENTARY_HESS_TEST:
817  case GA_NODE_ELEMENTARY_DIVERG_TEST:
818  is_elementary = true;
819  str << "Elementary_transformation(";
820  break;
821  case GA_NODE_SECONDARY_DOMAIN:
822  case GA_NODE_SECONDARY_DOMAIN_X:
823  case GA_NODE_SECONDARY_DOMAIN_NORMAL:
824  case GA_NODE_SECONDARY_DOMAIN_VAL:
825  case GA_NODE_SECONDARY_DOMAIN_GRAD:
826  case GA_NODE_SECONDARY_DOMAIN_HESS:
827  case GA_NODE_SECONDARY_DOMAIN_DIVERG:
828  case GA_NODE_SECONDARY_DOMAIN_VAL_TEST:
829  case GA_NODE_SECONDARY_DOMAIN_GRAD_TEST:
830  case GA_NODE_SECONDARY_DOMAIN_HESS_TEST:
831  case GA_NODE_SECONDARY_DOMAIN_DIVERG_TEST:
832  str << "Secondary_domain(";
833  is_secondary = true;
834  break;
835 
836  case GA_NODE_XFEM_PLUS:
837  case GA_NODE_XFEM_PLUS_VAL:
838  case GA_NODE_XFEM_PLUS_GRAD:
839  case GA_NODE_XFEM_PLUS_HESS:
840  case GA_NODE_XFEM_PLUS_DIVERG:
841  case GA_NODE_XFEM_PLUS_VAL_TEST:
842  case GA_NODE_XFEM_PLUS_GRAD_TEST:
843  case GA_NODE_XFEM_PLUS_HESS_TEST:
844  case GA_NODE_XFEM_PLUS_DIVERG_TEST:
845  is_xfem_plus = true;
846  str << "Xfem_plus(";
847  break;
848  case GA_NODE_XFEM_MINUS:
849  case GA_NODE_XFEM_MINUS_VAL:
850  case GA_NODE_XFEM_MINUS_GRAD:
851  case GA_NODE_XFEM_MINUS_HESS:
852  case GA_NODE_XFEM_MINUS_DIVERG:
853  case GA_NODE_XFEM_MINUS_VAL_TEST:
854  case GA_NODE_XFEM_MINUS_GRAD_TEST:
855  case GA_NODE_XFEM_MINUS_HESS_TEST:
856  case GA_NODE_XFEM_MINUS_DIVERG_TEST:
857  is_xfem_minus = true;
858  str << "Xfem_minus(";
859  break;
860  default:
861  break;
862  }
863 
864  switch(pnode->node_type) {
865  case GA_NODE_GRAD:
866  case GA_NODE_INTERPOLATE_GRAD:
867  case GA_NODE_ELEMENTARY_GRAD:
868  case GA_NODE_SECONDARY_DOMAIN_GRAD:
869  case GA_NODE_XFEM_PLUS_GRAD:
870  case GA_NODE_XFEM_MINUS_GRAD:
871  case GA_NODE_GRAD_TEST:
872  case GA_NODE_INTERPOLATE_GRAD_TEST:
873  case GA_NODE_ELEMENTARY_GRAD_TEST:
874  case GA_NODE_SECONDARY_DOMAIN_GRAD_TEST:
875  case GA_NODE_XFEM_PLUS_GRAD_TEST:
876  case GA_NODE_XFEM_MINUS_GRAD_TEST:
877  str << "Grad_";
878  break;
879  case GA_NODE_HESS:
880  case GA_NODE_INTERPOLATE_HESS:
881  case GA_NODE_ELEMENTARY_HESS:
882  case GA_NODE_SECONDARY_DOMAIN_HESS:
883  case GA_NODE_XFEM_PLUS_HESS:
884  case GA_NODE_XFEM_MINUS_HESS:
885  case GA_NODE_HESS_TEST:
886  case GA_NODE_INTERPOLATE_HESS_TEST:
887  case GA_NODE_ELEMENTARY_HESS_TEST:
888  case GA_NODE_SECONDARY_DOMAIN_HESS_TEST:
889  case GA_NODE_XFEM_PLUS_HESS_TEST:
890  case GA_NODE_XFEM_MINUS_HESS_TEST:
891  str << "Hess_";
892  break;
893  case GA_NODE_DIVERG:
894  case GA_NODE_INTERPOLATE_DIVERG:
895  case GA_NODE_SECONDARY_DOMAIN_DIVERG:
896  case GA_NODE_ELEMENTARY_DIVERG:
897  case GA_NODE_XFEM_PLUS_DIVERG:
898  case GA_NODE_XFEM_MINUS_DIVERG:
899  case GA_NODE_DIVERG_TEST:
900  case GA_NODE_INTERPOLATE_DIVERG_TEST:
901  case GA_NODE_ELEMENTARY_DIVERG_TEST:
902  case GA_NODE_SECONDARY_DOMAIN_DIVERG_TEST:
903  case GA_NODE_XFEM_PLUS_DIVERG_TEST:
904  case GA_NODE_XFEM_MINUS_DIVERG_TEST:
905  str << "Div_";
906  break;
907  default:
908  break;
909  }
910 
911  switch(pnode->node_type) {
912  case GA_NODE_OP:
913  {
914  bool par = false;
915  if (pnode->parent) {
916  if (pnode->parent->node_type == GA_NODE_OP &&
917  (ga_operator_priorities[pnode->op_type] >= 2 ||
918  ga_operator_priorities[pnode->op_type]
919  < ga_operator_priorities[pnode->parent->op_type]))
920  par = true;
921  if (pnode->parent->node_type == GA_NODE_PARAMS) par = true;
922  }
923 
924 
925  if (par) str << "(";
926  if (pnode->op_type == GA_UNARY_MINUS) {
927  GMM_ASSERT1(pnode->children.size() == 1, "Invalid tree");
928  str << "-"; ga_print_node(pnode->children[0], str);
929  } else if (pnode->op_type == GA_QUOTE) {
930  GMM_ASSERT1(pnode->children.size() == 1, "Invalid tree");
931  ga_print_node(pnode->children[0], str); str << "'";
932  } else if (pnode->op_type == GA_SYM) {
933  GMM_ASSERT1(pnode->children.size() == 1, "Invalid tree");
934  str << "Sym("; ga_print_node(pnode->children[0], str); str << ")";
935  } else if (pnode->op_type == GA_SKEW) {
936  GMM_ASSERT1(pnode->children.size() == 1, "Invalid tree");
937  str << "Skew("; ga_print_node(pnode->children[0], str); str << ")";
938  } else if (pnode->op_type == GA_TRACE) {
939  GMM_ASSERT1(pnode->children.size() == 1, "Invalid tree");
940  str << "Trace("; ga_print_node(pnode->children[0], str); str << ")";
941  } else if (pnode->op_type == GA_DEVIATOR) {
942  GMM_ASSERT1(pnode->children.size() == 1, "Invalid tree with "
943  << pnode->children.size() << " children instead of 1");
944  str << "Deviator("; ga_print_node(pnode->children[0], str); str<<")";
945  } else if (pnode->op_type == GA_PRINT) {
946  GMM_ASSERT1(pnode->children.size() == 1, "Invalid tree");
947  str << "Print("; ga_print_node(pnode->children[0], str); str << ")";
948  } else {
949  if (!par && pnode->op_type == GA_MULT &&
950  (pnode->children.size() == 1 ||
951  pnode->test_function_type == size_type(-1) ||
952  (pnode->children[0]->tensor_order() == 4 &&
953  pnode->children[1]->tensor_order() == 2)))
954  { par = true; str << "("; }
955  ga_print_node(pnode->children[0], str);
956  switch (pnode->op_type) {
957  case GA_PLUS: str << "+"; break;
958  case GA_MINUS: str << "-"; break;
959  case GA_MULT: str << "*"; break;
960  case GA_DIV: str << "/"; break;
961  case GA_COLON: str << ":"; break;
962  case GA_DOT: str << "."; break;
963  case GA_DOTMULT: str << ".*"; break;
964  case GA_DOTDIV: str << "./"; break;
965  case GA_TMULT: str << "@"; break;
966  default: GMM_ASSERT1(false, "Invalid or not taken into account "
967  "operation");
968  }
969  if (pnode->children.size() >= 2)
970  ga_print_node(pnode->children[1], str);
971  else
972  str << "(unknown second argument)";
973  }
974  if (par) str << ")";
975  }
976  break;
977 
978  case GA_NODE_X:
979  if (pnode->nbc1) str << "X(" << pnode->nbc1 << ")"; else str << "X";
980  break;
981  case GA_NODE_ELT_SIZE: str << "element_size"; break;
982  case GA_NODE_ELT_K: str << "element_K"; break;
983  case GA_NODE_ELT_B: str << "element_B"; break;
984  case GA_NODE_NORMAL: str << "Normal"; break;
985  case GA_NODE_INTERPOLATE_FILTER:
986  str << "Interpolate_filter(" << pnode->interpolate_name << ",";
987  ga_print_node(pnode->children[0], str);
988  if (pnode->children.size() == 2)
989  { str << ","; ga_print_node(pnode->children[1], str); }
990  else if (pnode->nbc1 != size_type(-1)) str << "," << pnode->nbc1;
991  str << ")";
992  break;
993  case GA_NODE_INTERPOLATE_X: case GA_NODE_SECONDARY_DOMAIN_X:
994  str << "X";
995  break;
996  case GA_NODE_INTERPOLATE_NORMAL: case GA_NODE_SECONDARY_DOMAIN_NORMAL:
997  str << "Normal";
998  break;
999  case GA_NODE_INTERPOLATE_ELT_K:
1000  str << "element_K";
1001  break;
1002  case GA_NODE_INTERPOLATE_ELT_B:
1003  str << "element_B";
1004  break;
1005  case GA_NODE_INTERPOLATE_DERIVATIVE:
1006  str << (pnode->test_function_type == 1 ? "Test_" : "Test2_")
1007  << "Interpolate_derivative(" << pnode->interpolate_name_der << ","
1008  << pnode->name;
1009  if (pnode->interpolate_name.size())
1010  str << "," << pnode->interpolate_name;
1011  str << ")";
1012  break;
1013  case GA_NODE_INTERPOLATE:
1014  case GA_NODE_ELEMENTARY:
1015  case GA_NODE_SECONDARY_DOMAIN:
1016  case GA_NODE_XFEM_PLUS:
1017  case GA_NODE_XFEM_MINUS:
1018  case GA_NODE_VAL:
1019  case GA_NODE_INTERPOLATE_VAL:
1020  case GA_NODE_ELEMENTARY_VAL:
1021  case GA_NODE_SECONDARY_DOMAIN_VAL:
1022  case GA_NODE_XFEM_PLUS_VAL:
1023  case GA_NODE_XFEM_MINUS_VAL:
1024  case GA_NODE_GRAD:
1025  case GA_NODE_INTERPOLATE_GRAD:
1026  case GA_NODE_SECONDARY_DOMAIN_GRAD:
1027  case GA_NODE_ELEMENTARY_GRAD:
1028  case GA_NODE_XFEM_PLUS_GRAD:
1029  case GA_NODE_XFEM_MINUS_GRAD:
1030  case GA_NODE_HESS:
1031  case GA_NODE_INTERPOLATE_HESS:
1032  case GA_NODE_SECONDARY_DOMAIN_HESS:
1033  case GA_NODE_ELEMENTARY_HESS:
1034  case GA_NODE_XFEM_PLUS_HESS:
1035  case GA_NODE_XFEM_MINUS_HESS:
1036  case GA_NODE_DIVERG:
1037  case GA_NODE_INTERPOLATE_DIVERG:
1038  case GA_NODE_ELEMENTARY_DIVERG:
1039  case GA_NODE_SECONDARY_DOMAIN_DIVERG:
1040  case GA_NODE_XFEM_PLUS_DIVERG:
1041  case GA_NODE_XFEM_MINUS_DIVERG:
1042  str << pnode->name;
1043  break;
1044  case GA_NODE_VAL_TEST:
1045  case GA_NODE_INTERPOLATE_VAL_TEST:
1046  case GA_NODE_ELEMENTARY_VAL_TEST:
1047  case GA_NODE_SECONDARY_DOMAIN_VAL_TEST:
1048  case GA_NODE_XFEM_PLUS_VAL_TEST:
1049  case GA_NODE_XFEM_MINUS_VAL_TEST:
1050  case GA_NODE_GRAD_TEST:
1051  case GA_NODE_INTERPOLATE_GRAD_TEST:
1052  case GA_NODE_ELEMENTARY_GRAD_TEST:
1053  case GA_NODE_SECONDARY_DOMAIN_GRAD_TEST:
1054  case GA_NODE_XFEM_PLUS_GRAD_TEST:
1055  case GA_NODE_XFEM_MINUS_GRAD_TEST:
1056  case GA_NODE_HESS_TEST:
1057  case GA_NODE_INTERPOLATE_HESS_TEST:
1058  case GA_NODE_ELEMENTARY_HESS_TEST:
1059  case GA_NODE_SECONDARY_DOMAIN_HESS_TEST:
1060  case GA_NODE_XFEM_PLUS_HESS_TEST:
1061  case GA_NODE_XFEM_MINUS_HESS_TEST:
1062  case GA_NODE_DIVERG_TEST:
1063  case GA_NODE_INTERPOLATE_DIVERG_TEST:
1064  case GA_NODE_ELEMENTARY_DIVERG_TEST:
1065  case GA_NODE_SECONDARY_DOMAIN_DIVERG_TEST:
1066  case GA_NODE_XFEM_PLUS_DIVERG_TEST:
1067  case GA_NODE_XFEM_MINUS_DIVERG_TEST:
1068  str << (pnode->test_function_type == 1 ? "Test_" : "Test2_")
1069  << pnode->name;
1070  break;
1071  case GA_NODE_SPEC_FUNC: str << pnode->name; break;
1072  case GA_NODE_OPERATOR:
1073  case GA_NODE_PREDEF_FUNC:
1074  if (pnode->der1) {
1075  str << "Derivative_" << pnode->der1 << "_";
1076  if (pnode->der2) str << pnode->der2 << "_";
1077  }
1078  str << pnode->name; break;
1079  case GA_NODE_ZERO:
1080  GMM_ASSERT1(pnode->test_function_type != size_type(-1),
1081  "Internal error");
1082  if (pnode->test_function_type) str << "(";
1083  ga_print_constant_tensor(pnode, str);
1084  if (pnode->name_test1.size()) {
1085  GMM_ASSERT1(pnode->qdim1 > 0, "Internal error");
1086  if (pnode->qdim1 == 1)
1087  str << "*Test_" << pnode->name_test1;
1088  else {
1089  str << "*(Reshape(Test_" << pnode->name_test1 << ","
1090  << pnode->qdim1<< ")(1))";
1091  }
1092  }
1093  if (pnode->name_test2.size()) {
1094  GMM_ASSERT1(pnode->qdim2 > 0, "Internal error");
1095  if (pnode->qdim2 == 1)
1096  str << "*Test2_" << pnode->name_test2;
1097  else {
1098  str << "*(Reshape(Test2_" << pnode->name_test2 << ","
1099  << pnode->qdim2<< ")(1))";
1100  }
1101  }
1102  if (pnode->test_function_type) str << ")";
1103  break;
1104 
1105  case GA_NODE_CONSTANT:
1106  ga_print_constant_tensor(pnode, str);
1107  break;
1108 
1109  case GA_NODE_ALLINDICES:
1110  str << ":";
1111  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1112  break;
1113 
1114  case GA_NODE_PARAMS:
1115  GMM_ASSERT1(pnode->children.size(), "Invalid tree");
1116  ga_print_node(pnode->children[0], str);
1117  str << "(";
1118  for (size_type i = 1; i < pnode->children.size(); ++i) {
1119  if (i > 1) str << ", ";
1120  ga_print_node(pnode->children[i], str);
1121  }
1122  str << ")";
1123  break;
1124 
1125  case GA_NODE_NAME:
1126  str << pnode->name;
1127  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1128  break;
1129 
1130  case GA_NODE_MACRO_PARAM:
1131  if (pnode->nbc2 == 1) str << "Grad_";
1132  if (pnode->nbc2 == 2) str << "Hess_";
1133  if (pnode->nbc2 == 3) str << "Div_";
1134  if (pnode->nbc3 == 1) str << "Test_";
1135  if (pnode->nbc3 == 2) str << "Test2_";
1136  str << "P" << pnode->nbc1;
1137  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1138  break;
1139 
1140  case GA_NODE_RESHAPE:
1141  str << "Reshape";
1142  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1143  break;
1144 
1145  case GA_NODE_CROSS_PRODUCT:
1146  str << "Cross_product";
1147  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1148  break;
1149 
1150  case GA_NODE_SWAP_IND:
1151  str << "Swap_indices";
1152  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1153  break;
1154 
1155  case GA_NODE_IND_MOVE_LAST:
1156  str << "Index_move_last";
1157  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1158  break;
1159 
1160  case GA_NODE_CONTRACT:
1161  str << "Contract";
1162  GMM_ASSERT1(pnode->children.size() == 0, "Invalid tree");
1163  break;
1164 
1165  case GA_NODE_C_MATRIX:
1166  GMM_ASSERT1(pnode->children.size(), "Invalid tree");
1167  GMM_ASSERT1(pnode->nbc1 == pnode->tensor_order(), "Invalid C_MATRIX");
1168  switch (pnode->tensor_order()) {
1169  case 0:
1170  ga_print_node(pnode->children[0], str);
1171  break;
1172 
1173  case 1:
1174  str << "[";
1175  for (size_type i = 0; i < pnode->tensor_proper_size(0); ++i) {
1176  if (i != 0) str << ",";
1177  ga_print_node(pnode->children[i], str);
1178  }
1179  str << "]";
1180  break;
1181 
1182  case 2: case 3: case 4:
1183  {
1184  size_type ii(0);
1185  size_type n0 = pnode->tensor_proper_size(0);
1186  size_type n1 = pnode->tensor_proper_size(1);
1187  size_type n2 = ((pnode->tensor_order() > 2) ?
1188  pnode->tensor_proper_size(2) : 1);
1189  size_type n3 = ((pnode->tensor_order() > 3) ?
1190  pnode->tensor_proper_size(3) : 1);
1191  if (n3 > 1) str << "[";
1192  for (size_type l = 0; l < n3; ++l) {
1193  if (l != 0) str << ",";
1194  if (n2 > 1) str << "[";
1195  for (size_type k = 0; k < n2; ++k) {
1196  if (k != 0) str << ",";
1197  if (n1 > 1) str << "[";
1198  for (size_type j = 0; j < n1; ++j) {
1199  if (j != 0) str << ",";
1200  if (n0 > 1) str << "[";
1201  for (size_type i = 0; i < n0; ++i) {
1202  if (i != 0) str << ",";
1203  ga_print_node(pnode->children[ii++], str);
1204  }
1205  if (n0 > 1) str << "]";
1206  }
1207  if (n1 > 1) str << "]";
1208  }
1209  if (n2 > 1) str << "]";
1210  }
1211  if (n3 > 1) str << "]";
1212  }
1213  break;
1214 
1215  case 5: case 6: case 7: case 8:
1216  str << "Reshape([";
1217  for (size_type i = 0; i < pnode->tensor_proper_size(); ++i) {
1218  if (i != 0) str << ";";
1219  ga_print_node(pnode->children[i], str);
1220  }
1221  str << "]";
1222  for (size_type i = 0; i < pnode->tensor_order(); ++i) {
1223  if (i != 0) str << ", ";
1224  str << pnode->tensor_proper_size(i);
1225  }
1226  str << ")";
1227  break;
1228 
1229  default: GMM_ASSERT1(false, "Invalid tensor dimension");
1230  }
1231  break;
1232 
1233  default:
1234  str << "Invalid or not taken into account node type "
1235  << pnode->node_type;
1236  break;
1237  }
1238 
1239  if (is_interpolate)
1240  str << "," << pnode->interpolate_name << ")";
1241  else if (is_elementary) {
1242  str << "," << pnode->elementary_name;
1243  if (pnode->name.compare(pnode->elementary_target) != 0)
1244  str << "," << pnode->elementary_target;
1245  str << ")";
1246  } else if (is_secondary)
1247  str << ")";
1248  else if (is_xfem_plus || is_xfem_minus)
1249  str << ")";
1250 
1251  str.precision(prec);
1252  }
1253 
1254  std::string ga_tree_to_string(const ga_tree &tree) {
1255  std::stringstream str;
1256  str.precision(16);
1257  if (tree.root) verify_tree(tree.root, 0);
1258  if (tree.root) ga_print_node(tree.root, str); else str << "0";
1259  return str.str();
1260  }
1261 
1262  size_type ga_parse_prefix_operator(std::string &name) {
1263  if (name.size() >= 5 && name.compare(0, 5, "Grad_") == 0)
1264  { name = name.substr(5); return 1; }
1265  else if (name.size() >= 5 && name.compare(0, 5, "Hess_") == 0)
1266  { name = name.substr(5); return 2; }
1267  else if (name.size() >= 4 && name.compare(0, 4, "Div_") == 0)
1268  { name = name.substr(4); return 3; }
1269  return 0;
1270  }
1271 
1272  size_type ga_parse_prefix_test(std::string &name) {
1273  if (name.size() >= 5 && name.compare(0, 5, "Test_") == 0)
1274  { name = name.substr(5); return 1; }
1275  else if (name.size() >= 6 && name.compare(0, 6, "Test2_") == 0)
1276  { name = name.substr(6); return 2; }
1277  return 0;
1278  }
1279 
1280  // 0 : ok
1281  // 1 : function or operator name or "X"
1282  // 2 : reserved prefix Grad, Hess, Div, Derivative_ Test and Test2
1283  // 3 : reserved prefix Dot and Previous
1284  int ga_check_name_validity(const std::string &name) {
1285  if (name.compare(0, 11, "Derivative_") == 0)
1286  return 2;
1287 
1288  const ga_predef_operator_tab &PREDEF_OPERATORS
1289  = dal::singleton<ga_predef_operator_tab>::instance(0);
1290  const ga_spec_function_tab &SPEC_FUNCTIONS
1291  = dal::singleton<ga_spec_function_tab>::instance(0);
1292  const ga_spec_op_tab &SPEC_OP
1293  = dal::singleton<ga_spec_op_tab>::instance(0);
1294  const ga_predef_function_tab &PREDEF_FUNCTIONS
1295  = dal::singleton<ga_predef_function_tab>::instance(0);
1296 
1297  if (SPEC_OP.find(name) != SPEC_OP.end())
1298  return 1;
1299 
1300  if (PREDEF_FUNCTIONS.find(name) != PREDEF_FUNCTIONS.end())
1301  return 1;
1302 
1303  if (SPEC_FUNCTIONS.find(name) != SPEC_FUNCTIONS.end())
1304  return 1;
1305 
1306  if (PREDEF_OPERATORS.tab.find(name) != PREDEF_OPERATORS.tab.end())
1307  return 1;
1308 
1309  if (name.size() >= 5 && name.compare(0, 5, "Grad_") == 0)
1310  return 2;
1311 
1312  if (name.size() >= 5 && name.compare(0, 5, "Hess_") == 0)
1313  return 2;
1314 
1315  if (name.size() >= 4 && name.compare(0, 4, "Div_") == 0)
1316  return 2;
1317 
1318  if (name.size() >= 6 && name.compare(0, 6, "Test2_") == 0)
1319  return 2;
1320 
1321  if (name.size() >= 5 && name.compare(0, 5, "Test_") == 0)
1322  return 2;
1323 
1324 // if (name.size() >= 4 && name.compare(0, 4, "Dot_") == 0)
1325 // return 3;
1326 // if (name.size() >= 5 && name.compare(0, 5, "Dot2_") == 0)
1327 // return 3;
1328 
1329 // if (name.size() >= 9 && name.compare(0, 9, "Previous_") == 0)
1330 // return 3;
1331 // if (name.size() >= 10 && name.compare(0, 10, "Previous2_") == 0)
1332 // return 3;
1333 // if (name.size() >= 12 && name.compare(0, 12, "Previous1_2_") == 0)
1334 // return 3;
1335 
1336  return 0;
1337  }
1338 
1339  //=========================================================================
1340  // Structure dealing with macros.
1341  //=========================================================================
1342 
1343  ga_macro::ga_macro() : ptree(new ga_tree), nbp(0) {}
1344  ga_macro::~ga_macro() { delete ptree; }
1345  ga_macro::ga_macro(const std::string &name, const ga_tree &t, size_type nbp_)
1346  : ptree(new ga_tree(t)), macro_name_(name), nbp(nbp_) {}
1347  ga_macro::ga_macro(const ga_macro &gam)
1348  : ptree(new ga_tree(gam.tree())), macro_name_(gam.name()),
1349  nbp(gam.nb_params()) {}
1350  ga_macro &ga_macro::operator =(const ga_macro &gam) {
1351  delete ptree; ptree = new ga_tree(gam.tree());
1352  macro_name_ = gam.name();
1353  nbp = gam.nb_params();
1354  return *this;
1355  }
1356 
1357  static void ga_replace_macro_params
1358  (ga_tree &tree, pga_tree_node pnode,
1359  const std::vector<pga_tree_node> &children) {
1360  if (!pnode) return;
1361  for (pga_tree_node &child : pnode->children)
1362  ga_replace_macro_params(tree, child, children);
1363 
1364  if (pnode->node_type == GA_NODE_MACRO_PARAM) {
1365  size_type po = pnode->nbc2;
1366  size_type pt = pnode->nbc3;
1367  GMM_ASSERT1(pnode->nbc1+1 < children.size(), "Internal error");
1368  pga_tree_node pchild = children[pnode->nbc1+1];
1369 
1370  if (po || pt || pnode->op_type != GA_NAME) {
1371  if (!(pchild->children.empty()) || pchild->node_type != GA_NODE_NAME)
1372  ga_throw_error(pchild->expr, pchild->pos, "Error in macro "
1373  "expansion. Only variable name are allowed for macro "
1374  "parameter preceded by Grad_ Hess_ Test_ or Test2_ "
1375  "prefixes.");
1376  switch(pnode->op_type) {
1377  case GA_NAME : pnode->node_type = GA_NODE_NAME; break;
1378  case GA_INTERPOLATE : pnode->node_type = GA_NODE_INTERPOLATE; break;
1379  case GA_INTERPOLATE_DERIVATIVE :
1380  pnode->node_type = GA_NODE_INTERPOLATE_DERIVATIVE; break;
1381  case GA_ELEMENTARY : pnode->node_type = GA_NODE_ELEMENTARY; break;
1382  case GA_SECONDARY_DOMAIN :
1383  pnode->node_type = GA_NODE_SECONDARY_DOMAIN; break;
1384  case GA_XFEM_PLUS : pnode->node_type = GA_NODE_XFEM_PLUS; break;
1385  case GA_XFEM_MINUS: pnode->node_type = GA_NODE_XFEM_MINUS; break;
1386  default:break;
1387  }
1388  pnode->name = pchild->name;
1389  if (pt == 1) pnode->name = "Test_" + pnode->name;
1390  if (pt == 2) pnode->name = "Test2_" + pnode->name;
1391  if (po == 1) pnode->name = "Grad_" + pnode->name;
1392  if (po == 2) pnode->name = "Hess_" + pnode->name;
1393  if (po == 3) pnode->name = "Div_" + pnode->name;
1394  } else {
1395  pga_tree_node pnode_old = pnode;
1396  pnode = nullptr;
1397  tree.copy_node(pchild, pnode);
1398  pnode->parent = pnode_old->parent;
1399  if (pnode_old->parent)
1400  pnode_old->parent->replace_child(pnode_old, pnode);
1401  else
1402  tree.root = pnode;
1403  GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
1404  delete pnode_old;
1405  }
1406  }
1407  }
1408 
1409  static void ga_expand_macro(ga_tree &tree, pga_tree_node pnode,
1410  const ga_macro_dictionary &macro_dict) {
1411  if (!pnode) return;
1412 
1413  if (pnode->node_type == GA_NODE_PARAMS) {
1414 
1415  for (size_type i = 1; i < pnode->children.size(); ++i)
1416  ga_expand_macro(tree, pnode->children[i], macro_dict);
1417 
1418  if (pnode->children[0]->node_type != GA_NODE_NAME) {
1419  ga_expand_macro(tree, pnode->children[0], macro_dict);
1420  } else {
1421 
1422  if (macro_dict.macro_exists(pnode->children[0]->name)) {
1423 
1424  const ga_macro &gam = macro_dict.get_macro(pnode->children[0]->name);
1425 
1426  if (gam.nb_params()==0) { // Macro without parameters
1427  pga_tree_node pnode_old = pnode->children[0];
1428  pnode->children[0] = nullptr;
1429  tree.copy_node(gam.tree().root, pnode->children[0]);
1430  pnode->children[0]->parent = pnode_old->parent;
1431  GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
1432  delete pnode_old;
1433  } else { // Macro with parameters
1434  if (gam.nb_params()+1 != pnode->children.size())
1435  ga_throw_error(pnode->expr, pnode->pos,
1436  "Bad number of parameters in the use of macro '"
1437  << gam.name() << "'. Expected " << gam.nb_params()
1438  << " found " << pnode->children.size()-1 << ".");
1439 
1440  pga_tree_node pnode_old = pnode;
1441  pnode = nullptr;
1442  tree.copy_node(gam.tree().root, pnode);
1443  pnode->parent = pnode_old->parent;
1444  if (pnode_old->parent)
1445  pnode_old->parent->replace_child(pnode_old, pnode);
1446  else
1447  tree.root = pnode;
1448  ga_replace_macro_params(tree, pnode, pnode_old->children);
1449  tree.clear_node_rec(pnode_old);
1450  }
1451  }
1452  }
1453 
1454  } else if (pnode->node_type == GA_NODE_NAME &&
1455  macro_dict.macro_exists(pnode->name)) {
1456  // Macro without parameters
1457  const ga_macro &gam = macro_dict.get_macro(pnode->name);
1458  if (gam.nb_params() != 0)
1459  ga_throw_error(pnode->expr, pnode->pos,
1460  "Bad number of parameters in the use of macro '"
1461  << gam.name() << "'. Expected " << gam.nb_params()
1462  << " none found.");
1463 
1464  pga_tree_node pnode_old = pnode;
1465  pnode = nullptr;
1466  tree.copy_node(gam.tree().root, pnode);
1467  pnode->parent = pnode_old->parent;
1468  if (pnode_old->parent)
1469  pnode_old->parent->replace_child(pnode_old, pnode);
1470  else
1471  tree.root = pnode;
1472  GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
1473  delete pnode_old;
1474  } else {
1475  for (pga_tree_node &child : pnode->children)
1476  ga_expand_macro(tree, child, macro_dict);
1477  }
1478  }
1479 
1480  static void ga_mark_macro_params_rec(const pga_tree_node pnode,
1481  const std::vector<std::string> &params) {
1482  if (!pnode) return;
1483  for (pga_tree_node &child : pnode->children)
1484  ga_mark_macro_params_rec(child, params);
1485 
1486  if (pnode->node_type == GA_NODE_NAME ||
1487  pnode->node_type == GA_NODE_INTERPOLATE ||
1488  pnode->node_type == GA_NODE_ELEMENTARY ||
1489  pnode->node_type == GA_NODE_SECONDARY_DOMAIN ||
1490  pnode->node_type == GA_NODE_XFEM_PLUS ||
1491  pnode->node_type == GA_NODE_XFEM_MINUS) {
1492  std::string name = pnode->name;
1493  size_type po = ga_parse_prefix_operator(name);
1494  size_type pt = ga_parse_prefix_test(name);
1495 
1496  for (size_type i = 0; i < params.size(); ++i)
1497  if (name.compare(params[i]) == 0) {
1498  pnode->name = name;
1499  switch(pnode->node_type) {
1500  case GA_NODE_NAME : pnode->op_type = GA_NAME; break;
1501  case GA_NODE_INTERPOLATE : pnode->op_type = GA_INTERPOLATE; break;
1502  case GA_NODE_INTERPOLATE_DERIVATIVE :
1503  pnode->op_type = GA_INTERPOLATE_DERIVATIVE; break;
1504  case GA_NODE_ELEMENTARY : pnode->op_type = GA_ELEMENTARY; break;
1505  case GA_NODE_SECONDARY_DOMAIN :
1506  pnode->op_type = GA_SECONDARY_DOMAIN; break;
1507  case GA_NODE_XFEM_PLUS : pnode->op_type = GA_XFEM_PLUS; break;
1508  case GA_NODE_XFEM_MINUS: pnode->op_type = GA_XFEM_MINUS; break;
1509  default:break;
1510  }
1511  pnode->node_type = GA_NODE_MACRO_PARAM;
1512  pnode->nbc1 = i; pnode->nbc2 = po; pnode->nbc3 = pt;
1513  }
1514  }
1515  }
1516 
1517  static void ga_mark_macro_params(ga_macro &gam,
1518  const std::vector<std::string> &params,
1519  const ga_macro_dictionary &macro_dict) {
1520  if (gam.tree().root) {
1521  ga_mark_macro_params_rec(gam.tree().root, params);
1522  ga_expand_macro(gam.tree(), gam.tree().root, macro_dict);
1523  }
1524  }
1525 
1526  bool ga_macro_dictionary::macro_exists(const std::string &name) const {
1527  if (macros.find(name) != macros.end()) return true;
1528  if (parent && parent->macro_exists(name)) return true;
1529  return false;
1530  }
1531 
1532  const ga_macro &
1533  ga_macro_dictionary::get_macro(const std::string &name) const {
1534  auto it = macros.find(name);
1535  if (it != macros.end()) return it->second;
1536  if (parent) return parent->get_macro(name);
1537  GMM_ASSERT1(false, "Undefined macro");
1538  }
1539 
1540  void ga_macro_dictionary::add_macro(const ga_macro &gam)
1541  { macros[gam.name()] = gam; }
1542 
1543  void ga_macro_dictionary::add_macro(const std::string &name,
1544  const std::string &expr)
1545  { ga_tree tree; ga_read_string_reg("Def "+name+":="+expr, tree, *this); }
1546 
1547  void ga_macro_dictionary::del_macro(const std::string &name) {
1548  auto it = macros.find(name);
1549  GMM_ASSERT1(it != macros.end(), "Undefined macro (at this level)");
1550  macros.erase(it);
1551  }
1552 
1553 
1554  //=========================================================================
1555  // Syntax analysis for the generic assembly language
1556  //=========================================================================
1557 
1558  // Read a term with an (implicit) pushdown automaton.
1559  static GA_TOKEN_TYPE ga_read_term(pstring expr, size_type &pos,
1560  ga_tree &tree,
1561  ga_macro_dictionary &macro_dict) {
1562  size_type token_pos, token_length;
1563  GA_TOKEN_TYPE t_type;
1564  int state = 1; // 1 = reading term, 2 = reading after term
1565 
1566  for (;;) {
1567 
1568  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1569 
1570  switch (state) {
1571 
1572  case 1:
1573  switch (t_type) {
1574  case GA_SCALAR:
1575  {
1576  char *endptr; const char *nptr = &((*expr)[token_pos]);
1577  scalar_type s_read = ::strtod(nptr, &endptr);
1578  if (endptr == nptr)
1579  ga_throw_error(expr, token_pos, "Bad numeric format.");
1580  tree.add_scalar(s_read, token_pos, expr);
1581  }
1582  state = 2; break;
1583 
1584  case GA_COLON:
1585  tree.add_allindices(token_pos, expr);
1586  state = 2; break;
1587 
1588  case GA_NAME:
1589  tree.add_name(&((*expr)[token_pos]), token_length, token_pos, expr);
1590  state = 2; break;
1591 
1592  case GA_MINUS: // unary -
1593  tree.add_op(GA_UNARY_MINUS, token_pos, expr);
1594  state = 1; break;
1595 
1596  case GA_PLUS: // unary +
1597  state = 1; break;
1598 
1599  case GA_SYM:
1600  tree.add_op(GA_SYM, token_pos, expr);
1601  state = 1; break;
1602 
1603  case GA_SKEW:
1604  tree.add_op(GA_SKEW, token_pos, expr);
1605  state = 1; break;
1606 
1607  case GA_TRACE:
1608  tree.add_op(GA_TRACE, token_pos, expr);
1609  state = 1; break;
1610 
1611  case GA_DEVIATOR:
1612  tree.add_op(GA_DEVIATOR, token_pos, expr);
1613  state = 1; break;
1614 
1615  case GA_DEF:
1616  {
1617  ga_macro gam;
1618  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1619  if (t_type != GA_NAME)
1620  ga_throw_error(expr, pos,
1621  "Macro definition should begin with macro name");
1622  gam.name() = std::string(&((*expr)[token_pos]), token_length);
1623  if (ga_check_name_validity(gam.name()))
1624  ga_throw_error(expr, pos-1, "Invalid macro name.")
1625  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1626  std::vector<std::string> params;
1627  if (t_type == GA_LPAR) {
1628  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1629  while (t_type == GA_NAME) {
1630  params.push_back(std::string(&((*expr)[token_pos]),
1631  token_length));
1632  if (ga_check_name_validity(params.back()))
1633  ga_throw_error(expr, pos-1, "Invalid macro parameter name.");
1634  for (size_type i = 0; i+1 < params.size(); ++i)
1635  if (params.back().compare(params[i]) == 0)
1636  ga_throw_error(expr, pos-1,
1637  "Invalid repeated macro parameter name.");
1638  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1639  if (t_type == GA_COMMA)
1640  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1641  }
1642  if (t_type != GA_RPAR)
1643  ga_throw_error(expr, pos-1,
1644  "Missing right parenthesis in macro definition.");
1645  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1646  }
1647  if (t_type != GA_COLON_EQ)
1648  ga_throw_error(expr, pos-1, "Missing := for macro definition.");
1649 
1650  t_type = ga_read_term(expr, pos, gam.tree(), macro_dict);
1651  if (gam.tree().root)
1652  ga_expand_macro(gam.tree(), gam.tree().root, macro_dict);
1653  gam.nb_params() = params.size();
1654  if (params.size())
1655  ga_mark_macro_params(gam, params, macro_dict);
1656  macro_dict.add_macro(gam);
1657 
1658  // cout << "macro \"" << gam.name() << "\" registered with "
1659  // << gam.nb_params() << " params := "
1660  // << ga_tree_to_string(gam.tree()) << endl;
1661 
1662  if (t_type == GA_END) return t_type;
1663  else if (t_type != GA_SEMICOLON)
1664  ga_throw_error(expr, pos-1,
1665  "Syntax error at the end of macro definition.");
1666  state = 1;
1667  }
1668  break;
1669 
1670  case GA_INTERPOLATE:
1671  {
1672  tree.add_scalar(scalar_type(0), token_pos, expr);
1673  tree.current_node->node_type = GA_NODE_INTERPOLATE;
1674  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1675  if (t_type != GA_LPAR)
1676  ga_throw_error(expr, pos-1, "Missing interpolate arguments.");
1677  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1678  if (t_type != GA_NAME)
1679  ga_throw_error(expr, pos,
1680  "First argument of Interpolate should be a "
1681  "variable, test function, X or Normal.");
1682  tree.current_node->name = std::string(&((*expr)[token_pos]),
1683  token_length);
1684 
1685  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1686  if (t_type != GA_COMMA)
1687  ga_throw_error(expr, pos, "Bad format for Interpolate "
1688  "arguments.");
1689  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1690  if (t_type != GA_NAME)
1691  ga_throw_error(expr, pos,
1692  "Second argument of Interpolate should be a "
1693  "transformation name.");
1694  tree.current_node->interpolate_name
1695  = std::string(&((*expr)[token_pos]), token_length);
1696  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1697  if (t_type != GA_RPAR)
1698  ga_throw_error(expr, pos-1, "Missing a parenthesis after "
1699  "interpolate arguments.");
1700  state = 2;
1701  }
1702  break;
1703 
1704  case GA_INTERPOLATE_DERIVATIVE:
1705  {
1706  tree.add_scalar(scalar_type(0), token_pos, expr);
1707  tree.current_node->node_type = GA_NODE_INTERPOLATE_DERIVATIVE;
1708  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1709  if (t_type != GA_LPAR)
1710  ga_throw_error(expr, pos-1,
1711  "Missing Interpolate_derivative arguments.");
1712  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1713  if (t_type != GA_NAME)
1714  ga_throw_error(expr, pos,
1715  "First argument of Interpolate should the "
1716  "interpolate transformation name ");
1717  tree.current_node->interpolate_name_der
1718  = std::string(&((*expr)[token_pos]), token_length);
1719  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1720  if (t_type != GA_COMMA)
1721  ga_throw_error(expr, pos, "Bad format for Interpolate_derivative "
1722  "arguments.");
1723  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1724  if (t_type != GA_NAME)
1725  ga_throw_error(expr, pos,
1726  "Second argument of Interpolate should be a "
1727  "variable name.");
1728  tree.current_node->name
1729  = std::string(&((*expr)[token_pos]), token_length);
1730 
1731  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1732  tree.current_node->interpolate_name = "";
1733  if (t_type == GA_COMMA) {
1734  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1735  if (t_type != GA_NAME)
1736  ga_throw_error(expr, pos,
1737  "Third argument of Interpolate should be a "
1738  "interpolate transformation name.");
1739  tree.current_node->interpolate_name
1740  = std::string(&((*expr)[token_pos]), token_length);
1741  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1742  }
1743  if (t_type != GA_RPAR)
1744  ga_throw_error(expr, pos-1, "Missing a parenthesis after "
1745  "Interpolate_derivative arguments.");
1746  state = 2;
1747  }
1748  break;
1749 
1750  case GA_ELEMENTARY:
1751  {
1752  tree.add_scalar(scalar_type(0), token_pos, expr);
1753  tree.current_node->node_type = GA_NODE_ELEMENTARY;
1754  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1755  if (t_type != GA_LPAR)
1756  ga_throw_error(expr, pos-1,
1757  "Missing Elementary_transformation arguments.");
1758  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1759  if (t_type != GA_NAME)
1760  ga_throw_error(expr, pos,
1761  "First argument of Elementary_transformation "
1762  "should be a variable or a test function.");
1763  tree.current_node->name = std::string(&((*expr)[token_pos]),
1764  token_length);
1765  tree.current_node->elementary_target = tree.current_node->name;
1766 
1767  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1768  if (t_type != GA_COMMA)
1769  ga_throw_error(expr, pos, "Bad format for "
1770  "Elementary_transformation arguments.");
1771  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1772  if (t_type != GA_NAME)
1773  ga_throw_error(expr, pos,
1774  "Second argument of Elementary_transformation "
1775  "should be a transformation name.");
1776  tree.current_node->elementary_name
1777  = std::string(&((*expr)[token_pos]), token_length);
1778  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1779 
1780  if (t_type == GA_COMMA) {
1781  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1782  if (t_type != GA_NAME)
1783  ga_throw_error(expr, pos,
1784  "Third argument of Elementary_transformation "
1785  "should be a variable or data name.");
1786 
1787  tree.current_node->elementary_target =
1788  std::string(&((*expr)[token_pos]), token_length);
1789  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1790  }
1791 
1792  if (t_type != GA_RPAR)
1793  ga_throw_error(expr, pos-1, "Missing a parenthesis after "
1794  "Elementary_transformation arguments.");
1795  state = 2;
1796  }
1797  break;
1798 
1799  case GA_SECONDARY_DOMAIN:
1800  {
1801  tree.add_scalar(scalar_type(0), token_pos, expr);
1802  tree.current_node->node_type = GA_NODE_SECONDARY_DOMAIN;
1803  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1804  if (t_type != GA_LPAR)
1805  ga_throw_error(expr, pos-1,"Missing Secondary_domain arguments.");
1806  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1807  if (t_type != GA_NAME)
1808  ga_throw_error(expr, pos,
1809  "First argument of Secondary_domain should be a "
1810  "variable, test function, X or Normal.");
1811  tree.current_node->name = std::string(&((*expr)[token_pos]),
1812  token_length);
1813  tree.current_node->interpolate_name = tree.secondary_domain;
1814  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1815  if (t_type != GA_RPAR)
1816  ga_throw_error(expr, pos-1, "Missing a parenthesis after "
1817  "Secondary_domain arguments.");
1818  state = 2;
1819  }
1820  break;
1821 
1822  case GA_XFEM_PLUS:
1823  {
1824  tree.add_scalar(scalar_type(0), token_pos, expr);
1825  tree.current_node->node_type = GA_NODE_XFEM_PLUS;
1826  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1827  if (t_type != GA_LPAR)
1828  ga_throw_error(expr, pos-1,
1829  "Missing Xfem_plus arguments.");
1830  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1831  if (t_type != GA_NAME)
1832  ga_throw_error(expr, pos,
1833  "The argument of Xfem_plus should be a "
1834  "variable or a test function.");
1835  tree.current_node->name = std::string(&((*expr)[token_pos]),
1836  token_length);
1837  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1838  if (t_type != GA_RPAR)
1839  ga_throw_error(expr, pos-1, "Missing a parenthesis after "
1840  "Xfem_plus argument.");
1841  state = 2;
1842  }
1843  break;
1844 
1845  case GA_XFEM_MINUS:
1846  {
1847  tree.add_scalar(scalar_type(0), token_pos, expr);
1848  tree.current_node->node_type = GA_NODE_XFEM_MINUS;
1849  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1850  if (t_type != GA_LPAR)
1851  ga_throw_error(expr, pos-1,
1852  "Missing Xfem_minus arguments.");
1853  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1854  if (t_type != GA_NAME)
1855  ga_throw_error(expr, pos,
1856  "The argument of Xfem_minus should be a "
1857  "variable or a test function.");
1858  tree.current_node->name = std::string(&((*expr)[token_pos]),
1859  token_length);
1860  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1861  if (t_type != GA_RPAR)
1862  ga_throw_error(expr, pos-1, "Missing a parenthesis after "
1863  "Xfem_minus argument.");
1864  state = 2;
1865  }
1866  break;
1867 
1868  case GA_INTERPOLATE_FILTER:
1869  {
1870  tree.add_scalar(scalar_type(0), token_pos, expr);
1871  tree.current_node->node_type = GA_NODE_INTERPOLATE_FILTER;
1872  tree.current_node->nbc1 = size_type(-1);
1873  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1874  if (t_type != GA_LPAR)
1875  ga_throw_error(expr, pos-1, "Missing interpolate arguments.");
1876  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1877  if (t_type != GA_NAME)
1878  ga_throw_error(expr, pos, "First argument of Interpolate_filter "
1879  "should be a transformation name.");
1880  tree.current_node->interpolate_name
1881  = std::string(&((*expr)[token_pos]), token_length);
1882  t_type = ga_get_token(*expr, pos, token_pos, token_length);
1883  if (t_type != GA_COMMA)
1884  ga_throw_error(expr, pos,
1885  "Bad format for Interpolate_filter arguments.");
1886  ga_tree sub_tree;
1887  t_type = ga_read_term(expr, pos, sub_tree, macro_dict);
1888  if (t_type != GA_RPAR && t_type != GA_COMMA)
1889  ga_throw_error(expr, pos-1,
1890  "Bad format for Interpolate_filter arguments.");
1891  tree.add_sub_tree(sub_tree);
1892  if (t_type == GA_COMMA) {
1893  ga_tree sub_tree2;
1894  t_type = ga_read_term(expr, pos, sub_tree2, macro_dict);
1895  tree.add_sub_tree(sub_tree2);
1896  }
1897  if (t_type != GA_RPAR)
1898  ga_throw_error(expr, pos-1, "Unbalanced parenthesis.");
1899  state = 2;
1900  }
1901  break;
1902 
1903  case GA_PRINT:
1904  tree.add_op(GA_PRINT, token_pos, expr);
1905  state = 1; break;
1906 
1907  case GA_LPAR: // Parenthesed expression
1908  {
1909  ga_tree sub_tree;
1910  GA_TOKEN_TYPE r_type;
1911  r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
1912  if (r_type != GA_RPAR)
1913  ga_throw_error(expr, pos-1, "Unbalanced parenthesis.");
1914  tree.add_sub_tree(sub_tree);
1915  state = 2;
1916  }
1917  break;
1918 
1919  case GA_LBRACKET: // Explicit vector/matrix or tensor
1920  {
1921  ga_tree sub_tree;
1922  GA_TOKEN_TYPE r_type;
1923  size_type nbc1(0), nbc2(0), nbc3(0), n1(0), n2(0), n3(0);
1924  size_type tensor_order(1);
1925  bool foundcomma(false), foundsemi(false);
1926 
1927  r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
1928  size_type nb_comp = 0;
1929  tree.add_matrix(token_pos, expr);
1930 
1931  if (sub_tree.root->node_type == GA_NODE_C_MATRIX) { // nested format
1932  bgeot::multi_index mii;
1933  do {
1934  if (nb_comp) {
1935  sub_tree.clear();
1936  r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
1937  }
1938  // in the nested format only "," and "]" are expected
1939  if (sub_tree.root->node_type != GA_NODE_C_MATRIX ||
1940  (r_type != GA_COMMA && r_type != GA_RBRACKET))
1941  ga_throw_error(expr, pos-1, "Bad explicit "
1942  "vector/matrix/tensor format.");
1943 
1944  // convert a row vector [a,b] to a column vector [a;b]
1945  if (sub_tree.root->marked &&
1946  sub_tree.root->tensor().sizes()[0] == 1 &&
1947  sub_tree.root->tensor().size() != 1) {
1948  bgeot::multi_index mi = sub_tree.root->tensor().sizes();
1949  for (size_type i = mi.size()-1; i > 0; i--)
1950  mi[i-1] = mi[i];
1951  mi.pop_back();
1952  sub_tree.root->tensor().adjust_sizes(mi);
1953  }
1954  if (!nb_comp)
1955  mii = sub_tree.root->tensor().sizes();
1956  else {
1957  const bgeot::multi_index &mi=sub_tree.root->tensor().sizes();
1958  bool cmp = true;
1959  if (mii.size() == mi.size()) {
1960  for (size_type i = 0; i < mi.size(); ++i)
1961  if (mi[i] != mii[i]) cmp = false;
1962  } else cmp = false;
1963  if (!cmp)
1964  ga_throw_error(expr, pos-1, "Bad explicit "
1965  "vector/matrix/tensor format.");
1966  }
1967  for (pga_tree_node &child : sub_tree.root->children) {
1968  child->parent = tree.current_node;
1969  tree.current_node->children.push_back(child);
1970  }
1971  sub_tree.root->children.resize(0);
1972  nb_comp++;
1973  } while (r_type != GA_RBRACKET);
1974  tree.current_node->marked = false;
1975  mii.push_back(nb_comp);
1976  tree.current_node->tensor().adjust_sizes(mii);
1977  } else { // non nested format
1978  do {
1979  if (nb_comp) {
1980  sub_tree.clear();
1981  r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
1982  }
1983  nb_comp++;
1984 
1985  tree.add_sub_tree(sub_tree);
1986 
1987  ++n1; ++n2; ++n3;
1988  if (tensor_order < 2) ++nbc1;
1989  if (tensor_order < 3) ++nbc2;
1990  if (tensor_order < 4) ++nbc3;
1991 
1992  if (r_type == GA_COMMA) {
1993  if (!foundcomma && tensor_order > 1)
1994  ga_throw_error(expr, pos-1, "Bad explicit "
1995  "vector/matrix/tensor format.");
1996  foundcomma = true;
1997  } else if (r_type == GA_SEMICOLON) {
1998  if (n1 != nbc1)
1999  ga_throw_error(expr, pos-1, "Bad explicit "
2000  "vector/matrix/tensor format.");
2001  n1 = 0;
2002  tensor_order = std::max(tensor_order, size_type(2));
2003  } else if (r_type == GA_DCOMMA) {
2004  if (n1 != nbc1 || n2 != nbc2)
2005  ga_throw_error(expr, pos-1, "Bad explicit "
2006  "vector/matrix/tensor format.");
2007  foundsemi = true;
2008  n2 = n1 = 0;
2009  tensor_order = std::max(tensor_order, size_type(3));
2010  } else if (r_type == GA_DSEMICOLON) {
2011  if (n1 != nbc1 || n2 != nbc2 || n3 != nbc3 ||
2012  tensor_order < 3)
2013  ga_throw_error(expr, pos-1, "Bad explicit "
2014  "vector/matrix/tensor format.");
2015  n3 = n2 = n1 = 0;
2016  tensor_order = std::max(tensor_order, size_type(4));
2017  } else if (r_type == GA_RBRACKET) {
2018  if (n1 != nbc1 || n2 != nbc2 || n3 != nbc3 ||
2019  tensor_order == 3)
2020  ga_throw_error(expr, pos-1, "Bad explicit "
2021  "vector/matrix/tensor format.");
2022  tree.current_node->nbc1 = nbc1;
2023  if (tensor_order == 4) {
2024  tree.current_node->nbc2 = nbc2/nbc1;
2025  tree.current_node->nbc3 = nbc3/nbc2;
2026  } else {
2027  tree.current_node->nbc2 = tree.current_node->nbc3 = 1;
2028  }
2029  } else {
2030  ga_throw_error(expr, pos-1, "The explicit "
2031  "vector/matrix/tensor components should be "
2032  "separated by ',', ';', ',,' and ';;' and "
2033  "be ended by ']'.");
2034  }
2035 
2036  } while (r_type != GA_RBRACKET);
2037  bgeot::multi_index mi;
2038  nbc1 = tree.current_node->nbc1;
2039  nbc2 = tree.current_node->nbc2;
2040  nbc3 = tree.current_node->nbc3;
2041 
2042  size_type nbl = tree.current_node->children.size()
2043  / (nbc2 * nbc1 * nbc3);
2044  switch(tensor_order) {
2045  case 1:
2046  // mi.push_back(1);
2047  mi.push_back(nbc1);
2048  break;
2049  case 2:
2050  mi.push_back(nbl);
2051  if (nbc1 > 1)
2052  mi.push_back(nbc1);
2053  break;
2054  case 3:
2055  mi.push_back(nbl);
2056  mi.push_back(nbc2);
2057  mi.push_back(nbc1);
2058  break;
2059  case 4:
2060  mi.push_back(nbl);
2061  mi.push_back(nbc3);
2062  mi.push_back(nbc2);
2063  mi.push_back(nbc1);
2064  break;
2065  default: GMM_ASSERT1(false, "Internal error");
2066  }
2067  tree.current_node->tensor().adjust_sizes(mi);
2068  std::vector<pga_tree_node> children = tree.current_node->children;
2069  auto it = tree.current_node->children.begin();
2070  for (size_type i = 0; i < nbc1; ++i)
2071  for (size_type j = 0; j < nbc2; ++j)
2072  for (size_type k = 0; k < nbc3; ++k)
2073  for (size_type l = 0; l < nbl; ++l, ++it)
2074  *it = children[i+nbc1*(j+nbc2*(k+nbc3*l))];
2075  tree.current_node->marked = true;
2076  }
2077  }
2078  tree.current_node->nbc1 = tree.current_node->tensor().sizes().size();
2079  state = 2;
2080  break;
2081 
2082  default:
2083  ga_throw_error(expr, token_pos, "Unexpected token.");
2084  }
2085  break;
2086 
2087  case 2:
2088  switch (t_type) {
2089  case GA_PLUS: case GA_MINUS: case GA_MULT: case GA_DIV:
2090  case GA_COLON: case GA_DOT: case GA_DOTMULT: case GA_DOTDIV:
2091  case GA_TMULT:
2092  tree.add_op(t_type, token_pos, expr);
2093  state = 1;
2094  break;
2095  case GA_QUOTE:
2096  tree.add_op(t_type, token_pos, expr);
2097  state = 2;
2098  break;
2099  case GA_END: case GA_RPAR: case GA_COMMA: case GA_DCOMMA:
2100  case GA_RBRACKET: case GA_SEMICOLON: case GA_DSEMICOLON:
2101  return t_type;
2102  case GA_LPAR: // Parameter list
2103  {
2104  ga_tree sub_tree;
2105  GA_TOKEN_TYPE r_type;
2106  tree.add_params(token_pos, expr);
2107  do {
2108  r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
2109  if (r_type != GA_RPAR && r_type != GA_COMMA)
2110  ga_throw_error(expr, pos-((r_type != GA_END)?1:0),
2111  "Parameters should be separated "
2112  "by ',' and parameter list ended by ')'.");
2113  tree.add_sub_tree(sub_tree);
2114  } while (r_type != GA_RPAR);
2115  state = 2;
2116  }
2117  break;
2118 
2119  default:
2120  ga_throw_error(expr, token_pos, "Unexpected token.");
2121  }
2122  break;
2123  }
2124  }
2125 
2126  return GA_INVALID;
2127  }
2128 
2129  // Syntax analysis of a string. Conversion to a tree. register the macros.
2130  void ga_read_string_reg(const std::string &expr, ga_tree &tree,
2131  ga_macro_dictionary &macro_dict) {
2132  size_type pos = 0, token_pos, token_length;
2133  tree.clear();
2134  GA_TOKEN_TYPE t = ga_get_token(expr, pos, token_pos, token_length);
2135  if (t == GA_END) return;
2136  pos = 0;
2137  pstring nexpr(new std::string(expr));
2138 
2139  t = ga_read_term(nexpr, pos, tree, macro_dict);
2140  if (tree.root) ga_expand_macro(tree, tree.root, macro_dict);
2141 
2142  switch (t) {
2143  case GA_RPAR:
2144  ga_throw_error(nexpr, pos-1, "Unbalanced parenthesis.");
2145  break;
2146  case GA_RBRACKET:
2147  ga_throw_error(nexpr, pos-1, "Unbalanced bracket.");
2148  break;
2149  case GA_END:
2150  break;
2151  default:
2152  ga_throw_error(nexpr, pos-1, "Unexpected token.");
2153  break;
2154  }
2155  }
2156 
2157  // Syntax analysis of a string. Conversion to a tree.
2158  // Do not register the macros (but expand them).
2159  void ga_read_string(const std::string &expr, ga_tree &tree,
2160  const ga_macro_dictionary &macro_dict) {
2161  ga_macro_dictionary macro_dict_loc(true, macro_dict);
2162  ga_read_string_reg(expr, tree, macro_dict_loc);
2163  }
2164 
2165  // Small tool to make basic substitutions into an assembly string
2166  // Should be replaced by macros now.
2167  std::string ga_substitute(const std::string &expr,
2168  const std::map<std::string, std::string> &dict) {
2169  if (dict.size()) {
2170  size_type pos = 0, token_pos, token_length;
2171  std::stringstream exprs;
2172 
2173  while (true) {
2174  GA_TOKEN_TYPE t_type = ga_get_token(expr, pos, token_pos, token_length);
2175  if (t_type == GA_END) return exprs.str();
2176  std::string name(&(expr[token_pos]), token_length);
2177  if (t_type == GA_NAME && dict.find(name) != dict.end())
2178  exprs << dict.at(name); else exprs << name;
2179  }
2180  }
2181  return expr;
2182  }
2183 
2184 } /* end of namespace */
dal::singleton::instance
static T & instance()
Instance from the current thread.
Definition: dal_singleton.h:165
getfem_generic_assembly_tree.h
Compilation and execution operations.
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
getfem
GEneric Tool for Finite Element Methods.
Definition: getfem_accumulated_distro.h:47