[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Getfem-commits] (no subject)
|
From: |
Yves Renard |
|
Subject: |
[Getfem-commits] (no subject) |
|
Date: |
Tue, 27 Aug 2019 14:37:03 -0400 (EDT) |
branch: master
commit bb28ea4115758c7de43f1b80790a1851de41d22f
Author: Yves Renard <address@hidden>
Date: Thu Aug 15 13:38:17 2019 +0200
Allowing elementary transformations between two different fems
---
doc/sphinx/source/userdoc/gasm_high.rst | 23 +-
.../source/userdoc/model_generic_assembly.rst | 2 +-
interface/tests/matlab/Makefile.am | 1 +
interface/tests/matlab/check_all.m | 8 +
interface/tests/matlab/check_mitc.m | 121 +++
src/getfem/getfem_generic_assembly.h | 4 +-
.../getfem_generic_assembly_compile_and_exec.h | 6 +-
src/getfem/getfem_generic_assembly_tree.h | 4 +-
src/getfem_generic_assembly_compile_and_exec.cc | 810 +++++++++++----------
src/getfem_generic_assembly_semantic.cc | 21 +-
src/getfem_generic_assembly_tree.cc | 25 +-
src/getfem_linearized_plates.cc | 16 +-
12 files changed, 626 insertions(+), 415 deletions(-)
diff --git a/doc/sphinx/source/userdoc/gasm_high.rst
b/doc/sphinx/source/userdoc/gasm_high.rst
index 24c52d9..526115f 100644
--- a/doc/sphinx/source/userdoc/gasm_high.rst
+++ b/doc/sphinx/source/userdoc/gasm_high.rst
@@ -70,7 +70,7 @@ A specific weak form language has been developed to describe
the weak formulatio
- ``Interpolate(variable, transformation)``: Powerful operation which allows
to interpolate the variables, or test functions either on the same mesh on
other elements or on another mesh. ``transformation`` is an object stored by
the workspace or model object which describes the map from the current point to
the point where to perform the interpolation. This functionality can be used
for instance to prescribe periodic conditions or to compute mortar matrices for
two finite element space [...]
- - ``Elementary_transformation(variable, transformation)``: Allow a linear
tranformation defined at the element level (i.e. not possible to define at the
gauss point level). This feature has been added mostly for defining a reduction
for plate elements (projection onto low-level vector element such as rotated
RT0). ``transformation`` is an object stored by the workspace or model object
which describes the trasformation for a particular element.
+ - ``Elementary_transformation(variable, transformation, dest)``: Allow a
linear tranformation defined at the element level (i.e. not possible to define
at the gauss point level). This feature has been added mostly for defining a
reduction for plate elements (projection onto low-level vector element such as
rotated RT0). ``transformation`` is an object stored by the workspace or model
object which describes the trasformation for a particular element. ``dest`` is
an optional argument ref [...]
- Possibility of integration on the direct product of two-domains for double
integral computation or coupling of two variables with a Kernel / convolution /
exchange integral. This allows terms like
:math:`\displaystyle\int_{\Omega_1}\int_{\Omega_2}k(x,y)u(x)v(y)dydx` with
:math:`\Omega_1` and :math:`\Omega_2` two domains, different or not, having
their own meshes, integration methods and with :math:`u` a variable defined on
:math:`\Omega_1` and :math:`v` a variable defined on :math:`\ [...]
@@ -918,16 +918,19 @@ the method::
where ``pelementary_transformation`` is a pointer to an object deriving from
``virtual_elementary_transformation``. Once it is added to the model/workspace,
it is possible to use the following expressions in the weak form language::
- Elementary_transformation(u, transname)
- Elementary_transformation(Grad_u, transname)
- Elementary_transformation(Div_u, transname)
- Elementary_transformation(Hess_u, transname)
- Elementary_transformation(Test_u, transname)
- Elementary_transformation(Grad_Test_u, transname)
- Elementary_transformation(Div_Test_u, transname)
- Elementary_transformation(Hess_Test_u, transname)
+ Elementary_transformation(u, transname[, dest])
+ Elementary_transformation(Grad_u, transname[, dest])
+ Elementary_transformation(Div_u, transname[, dest])
+ Elementary_transformation(Hess_u, transname[, dest])
+ Elementary_transformation(Test_u, transname[, dest])
+ Elementary_transformation(Grad_Test_u, transname[, dest])
+ Elementary_transformation(Div_Test_u, transname[, dest])
+ Elementary_transformation(Hess_Test_u, transname[, dest])
-where ``u`` is one of the FEM variables of the model/workspace. For the
moment, the only available elementary transformation is the the one for the
projection on rotated RT0 element for two-dimensional elements which can be
added thanks to the function (defined in
:file:`src/getfem/getfem_linearized_plates.h`)::
+where ``u`` is one of the FEM variables of the model/workspace, and ``dest``
is an optional parameter which should be a variable or data name of the model
and will correspond to the target fem of the transformation. If omitted, by
default, the transformation is from the fem of the first variable to itself.
+
+
+For the moment, the only available elementary transformation is the the one
for the projection on rotated RT0 element for two-dimensional elements which
can be added thanks to the function (defined in
:file:`src/getfem/getfem_linearized_plates.h`)::
add_2D_rotated_RT0_projection(model, transname)
diff --git a/doc/sphinx/source/userdoc/model_generic_assembly.rst
b/doc/sphinx/source/userdoc/model_generic_assembly.rst
index 64af4fe..82db0a4 100644
--- a/doc/sphinx/source/userdoc/model_generic_assembly.rst
+++ b/doc/sphinx/source/userdoc/model_generic_assembly.rst
@@ -1,4 +1,4 @@
-.. $Id: model_generic_assembly.rst 3655 2010-07-17 20:42:08Z renard $
+.. $Id: model_generic_assembly.rst 3655 2010-07-19 20:42:08Z renard $
.. include:: ../replaces.txt
diff --git a/interface/tests/matlab/Makefile.am
b/interface/tests/matlab/Makefile.am
index 442c69a..b25a141 100644
--- a/interface/tests/matlab/Makefile.am
+++ b/interface/tests/matlab/Makefile.am
@@ -37,6 +37,7 @@ EXTRA_DIST= \
check_plot.m \
check_slices.m \
check_spmat.m \
+ check_mitc.m \
check_workspace.m \
check_plasticity.m \
demo_elasticity.m \
diff --git a/interface/tests/matlab/check_all.m
b/interface/tests/matlab/check_all.m
index 17764aa..cd80a8c 100644
--- a/interface/tests/matlab/check_all.m
+++ b/interface/tests/matlab/check_all.m
@@ -96,6 +96,14 @@ catch
errcnt=errcnt+1; disp(['== ' t ' : FAILURE']);
end;
+t = 'check_mitc [check mitc4 element and elementary transformations] ';
+try
+ check_mitc;
+ disp(['== ' t ': SUCCESS']);
+catch
+ errcnt=errcnt+1; disp(['== ' t ' : FAILURE']);
+end;
+
t = 'demo_laplacian [model use for solving a Poisson problem] ';
try
automatic_var654 = 1;
diff --git a/interface/tests/matlab/check_mitc.m
b/interface/tests/matlab/check_mitc.m
new file mode 100644
index 0000000..f03c4a7
--- /dev/null
+++ b/interface/tests/matlab/check_mitc.m
@@ -0,0 +1,121 @@
+% Copyright (C) 2015-2017 FABRE Mathieu, SECK Mamadou, DALLERIT Valentin,
+%
+% This file is a part of GetFEM++
+%
+% GetFEM++ is free software; you can redistribute it and/or modify it
+% under the terms of the GNU Lesser General Public License as published
+% by the Free Software Foundation; either version 3 of the License, or
+% (at your option) any later version along with the GCC Runtime Library
+% Exception either version 3.1 or (at your option) any later version.
+% This program is distributed in the hope that it will be useful, but
+% WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
+% or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
+% License and GCC Runtime Library Exception for more details.
+% You should have received a copy of the GNU Lesser General Public License
+% along with this program; if not, write to the Free Software Foundation,
+% Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
+
+% Check on a simple supported Mindlin-Reissner plate
+
+function check_mitc(iverbose,idebug)
+ global gverbose;
+ global gdebug;
+ if (nargin >= 1),
+ gverbose = iverbose;
+ if (nargin == 2),
+ gdebug = idebug;
+ else
+ gdebug = 0;
+ end;
+ else
+ gverbose = 0; gdebug = 0;
+ end;
+
+
+
+
+Emodulus = 1; % Young Modulus
+nu = 0.5; % Poisson Coefficient
+epsilon = 0.001; % Plate thickness
+kappa = 5/6; % Shear correction factor
+f = -5*epsilon^3; % Prescribed force on the top of the plate
+
+variant = 2; % 0 : not reduced, 1 : with reduced integration,
+ % 2 : MITC reduction
+quadrangles = true; % Locking free only on quadrangle for the moment
+K = 1; % Degree of the finite element method
+with_Mindlin_brick = true; % Uses the Reissner-Mindlin predefined brick or not
+dirichlet_version = 1; % 0 = simplification, 1 = with multipliers
+ % 2 = penalization
+r = 1E8; % Penalization parameter.
+
+% trace on;
+gf_workspace('clear all');
+NX = 80;
+if (quadrangles)
+ m = gf_mesh('cartesian',[0:1/NX:1],[0:1/NX:1]);
+else
+
m=gf_mesh('import','structured',sprintf('GT="GT_PK(2,1)";SIZES=[1,1];NOISED=0;NSUBDIV=[%d,%d];',
NX, NX));
+end
+
+% Create a mesh_fem for a 2D vector field
+mftheta = gf_mesh_fem(m,2);
+mfu = gf_mesh_fem(m,1);
+if (quadrangles)
+ gf_mesh_fem_set(mftheta,'fem',gf_fem(sprintf('FEM_QK(2,%d)', K)));
+ gf_mesh_fem_set(mfu,'fem',gf_fem(sprintf('FEM_QK(2,%d)', K)));
+ mim = gf_mesh_im(m, gf_integ('IM_GAUSS_PARALLELEPIPED(2,6)'));
+ mim_reduced = gf_mesh_im(m, gf_integ('IM_GAUSS_PARALLELEPIPED(2,1)'));
+else
+ gf_mesh_fem_set(mftheta,'fem',gf_fem(sprintf('FEM_PK(2,%d)', K)));
+ gf_mesh_fem_set(mfu,'fem',gf_fem(sprintf('FEM_PK(2,%d)', K)));
+ mim = gf_mesh_im(m, gf_integ('IM_TRIANGLE(6)'));
+ mim_reduced = gf_mesh_im(m, gf_integ('IM_TRIANGLE(1)'));
+end
+
+% Detect the border of the mesh and assign it the boundary number 1
+border = gf_mesh_get(m,'outer faces');
+gf_mesh_set(m, 'boundary', 1, border);
+
+% Build the model
+md=gf_model('real');
+gf_model_set(md, 'add fem variable', 'u', mfu);
+gf_model_set(md, 'add fem variable', 'theta', mftheta);
+gf_model_set(md, 'add initialized data', 'E', Emodulus);
+gf_model_set(md, 'add initialized data', 'nu', nu);
+gf_model_set(md, 'add initialized data', 'epsilon', epsilon);
+gf_model_set(md, 'add initialized data', 'kappa', kappa);
+
+
+if (with_Mindlin_brick)
+ gf_model_set(md, 'add Mindlin Reissner plate brick', mim, mim_reduced, 'u',
'theta', 'E', 'nu', 'epsilon', 'kappa', variant);
+else
+ gf_model_set(md, 'add elementary rotated RT0 projection', 'RT0_projection');
+ gf_model_set(md, 'add linear term', mim,
'(E*epsilon*epsilon*epsilon*(1-nu)/(48 * (1 - nu*nu))) *
((Grad_theta+Grad_theta''):(Grad_Test_theta+Grad_Test_theta''))');
+ gf_model_set(md, 'add linear term', mim, '(E*epsilon*epsilon*epsilon*nu/(12
* (1 - nu*nu))) * (Div_theta*Div_Test_theta)');
+ if (variant == 0)
+ gf_model_set(md, 'add linear term', mim, '(E*kappa*epsilon/(1 + nu)) *
((Grad_u + theta).Grad_Test_u) + (E*kappa*epsilon/(1 + nu)) * ((Grad_u +
theta).Test_theta)');
+ elseif (variant == 1)
+ gf_model_set(md, 'add linear term', mim_reduced, '(E*kappa*epsilon/(1 +
nu)) * ((Grad_u + theta).Grad_Test_u) + (E*kappa*epsilon/(1 + nu)) * ((Grad_u +
theta).Test_theta)');
+ else
+ gf_model_set(md, 'add linear term', mim, '(E*kappa*epsilon/(1 + nu)) *
((Grad_u + Elementary_transformation(theta,RT0_projection)).Grad_Test_u) +
(E*kappa*epsilon/(1 + nu)) * ((Grad_u + Elementary_transformation(theta,
RT0_projection)).(Elementary_transformation(Test_theta, RT0_projection)))');
+ end
+end
+
+gf_model_set(md, 'add initialized data', 'VolumicData', f);
+
+gf_model_set(md, 'add source term brick', mim, 'u', 'VolumicData');
+gf_model_set(md, 'add initialized data', 'DirichletData', 0);
+switch (dirichlet_version)
+ case 0,
+ gf_model_set(md, 'add Dirichlet condition with simplification', 'u', 1,
'DirichletData');
+ case 1,
+ gf_model_set(md, 'add Dirichlet condition with multipliers', mim, 'u',
mfu, 1, 'DirichletData');
+ case 2,
+ gf_model_set(md, 'add Dirichlet condition with penalization', mim, 'u', r,
1, 'DirichletData');
+end
+gf_model_get(md, 'solve');
+U = gf_model_get(md, 'variable', 'u');
+
+gfassert('abs(min(U) + 0.1828) < 0.001');
+
diff --git a/src/getfem/getfem_generic_assembly.h
b/src/getfem/getfem_generic_assembly.h
index 3a05caa..68ea575 100644
--- a/src/getfem/getfem_generic_assembly.h
+++ b/src/getfem/getfem_generic_assembly.h
@@ -124,8 +124,8 @@ namespace getfem {
public:
- virtual void give_transformation(const mesh_fem &mf, size_type cv,
- base_matrix &M) const = 0;
+ virtual void give_transformation(const mesh_fem &mf1, const mesh_fem &mf2,
+ size_type cv, base_matrix &M) const = 0;
virtual ~virtual_elementary_transformation() {}
};
diff --git a/src/getfem/getfem_generic_assembly_compile_and_exec.h
b/src/getfem/getfem_generic_assembly_compile_and_exec.h
index 5c7ef03..ba49d6b 100644
--- a/src/getfem/getfem_generic_assembly_compile_and_exec.h
+++ b/src/getfem/getfem_generic_assembly_compile_and_exec.h
@@ -1,6 +1,6 @@
/*===========================================================================
- Copyright (C) 2013-2018 Yves Renard
+ Copyright (C) 2013-2019 Yves Renard
This file is a part of GetFEM++
@@ -146,7 +146,6 @@ namespace getfem {
struct elementary_trans_info {
base_matrix M;
- const mesh_fem *mf;
size_type icv;
};
@@ -178,7 +177,8 @@ namespace getfem {
std::map<std::string, std::set<std::string>> transformations;
std::set<std::string> transformations_der;
std::map<std::string, interpolate_info> interpolate_infos;
- std::map<std::string, elementary_trans_info> elementary_trans_infos;
+ std::map<std::tuple<std::string, const mesh_fem *, const mesh_fem *>,
+ elementary_trans_info> elementary_trans_infos;
secondary_domain_info secondary_domain_infos;
std::vector<pga_instruction>
diff --git a/src/getfem/getfem_generic_assembly_tree.h
b/src/getfem/getfem_generic_assembly_tree.h
index 6a9ae7f..7715299 100644
--- a/src/getfem/getfem_generic_assembly_tree.h
+++ b/src/getfem/getfem_generic_assembly_tree.h
@@ -1,6 +1,6 @@
/*===========================================================================
- Copyright (C) 2013-2018 Yves Renard
+ Copyright (C) 2013-2019 Yves Renard
This file is a part of GetFEM++
@@ -327,6 +327,8 @@ namespace getfem {
// name of transformation
std::string elementary_name; // For Elementary_transformation :
// name of transformation
+ std::string elementary_target;// For Elementary_transformation :
+ // target variable (for its mesh_fem)
size_type der1, der2; // For functions and nonlinear operators,
// optional derivative or second derivative.
bool symmetric_op;
diff --git a/src/getfem_generic_assembly_compile_and_exec.cc
b/src/getfem_generic_assembly_compile_and_exec.cc
index cbf0637..3614cd9 100644
--- a/src/getfem_generic_assembly_compile_and_exec.cc
+++ b/src/getfem_generic_assembly_compile_and_exec.cc
@@ -1,6 +1,6 @@
/*===========================================================================
- Copyright (C) 2013-2018 Yves Renard
+ Copyright (C) 2013-2019 Yves Renard
This file is a part of GetFEM++
@@ -1105,111 +1105,110 @@ namespace getfem {
: ga_instruction_copy_val_base(tt, Z_, q) {}
};
- struct ga_instruction_elementary_transformation {
+ struct ga_instruction_elementary_trans {
const base_vector &coeff_in;
base_vector coeff_out;
pelementary_transformation elemtrans;
- const mesh_fem &mf;
+ const mesh_fem &mf1, &mf2;
const fem_interpolation_context &ctx;
base_matrix &M;
- const mesh_fem **mf_M;
size_type &icv;
- void do_transformation() {
- size_type nn = gmm::vect_size(coeff_in);
- if (M.size() == 0 || icv != ctx.convex_num() || &mf != *mf_M) {
- M.base_resize(nn, nn);
- *mf_M = &mf; icv = ctx.convex_num();
- elemtrans->give_transformation(mf, icv, M);
+ void do_transformation(size_type n, size_type m) {
+ if (icv != ctx.convex_num() || M.size() == 0) {
+ M.base_resize(m, n);
+ icv = ctx.convex_num();
+ elemtrans->give_transformation(mf1, mf2, icv, M);
}
- coeff_out.resize(nn);
+ coeff_out.resize(gmm::mat_nrows(M));
gmm::mult(M, coeff_in, coeff_out); // remember: coeff == coeff_out
}
- ga_instruction_elementary_transformation
+ ga_instruction_elementary_trans
(const base_vector &co, pelementary_transformation e,
- const mesh_fem &mf_, const fem_interpolation_context &ctx_,
- base_matrix &M_, const mesh_fem **mf_M_, size_type &icv_)
- : coeff_in(co), elemtrans(e), mf(mf_), ctx(ctx_),
- M(M_), mf_M(mf_M_), icv(icv_) {}
- ~ga_instruction_elementary_transformation() {};
+ const mesh_fem &mf1_, const mesh_fem &mf2_,
+ const fem_interpolation_context &ctx_, base_matrix &M_,
+ size_type &icv_)
+ : coeff_in(co), elemtrans(e), mf1(mf1_), mf2(mf2_), ctx(ctx_),
+ M(M_), icv(icv_) {}
+ ~ga_instruction_elementary_trans() {};
};
- struct ga_instruction_elementary_transformation_val
- : public ga_instruction_val, ga_instruction_elementary_transformation {
+ struct ga_instruction_elementary_trans_val
+ : public ga_instruction_val, ga_instruction_elementary_trans {
// Z(ndof,target_dim), coeff_in(Qmult,ndof) --> t(target_dim*Qmult)
virtual int exec() {
GA_DEBUG_INFO("Instruction: variable value with elementary "
"transformation");
- do_transformation();
+ size_type ndof = Z.sizes()[0];
+ size_type Qmult = qdim / Z.sizes()[1];
+ do_transformation(ndof*Qmult, t.sizes()[0]);
return ga_instruction_val::exec();
}
- ga_instruction_elementary_transformation_val
+ ga_instruction_elementary_trans_val
(base_tensor &tt, const base_tensor &Z_, const base_vector &co, size_type
q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_val(tt, Z_, coeff_out, q),
- ga_instruction_elementary_transformation(co, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans(co, e, mf1_, mf2_, ctx_, M_, icv_) {}
};
- struct ga_instruction_elementary_transformation_grad
- : public ga_instruction_grad, ga_instruction_elementary_transformation {
+ struct ga_instruction_elementary_trans_grad
+ : public ga_instruction_grad, ga_instruction_elementary_trans {
// Z(ndof,target_dim,N), coeff_in(Qmult,ndof) --> t(target_dim*Qmult,N)
virtual int exec() {
GA_DEBUG_INFO("Instruction: gradient with elementary transformation");
- do_transformation();
+ size_type ndof = Z.sizes()[0];
+ size_type Qmult = qdim / Z.sizes()[1];
+ do_transformation(ndof*Qmult, t.sizes()[0]);
return ga_instruction_grad::exec();
}
- ga_instruction_elementary_transformation_grad
+ ga_instruction_elementary_trans_grad
(base_tensor &tt, const base_tensor &Z_, const base_vector &co, size_type
q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_grad(tt, Z_, coeff_out, q),
- ga_instruction_elementary_transformation(co, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans(co, e, mf1_, mf2_, ctx_, M_, icv_) {}
};
- struct ga_instruction_elementary_transformation_hess
- : public ga_instruction_hess, ga_instruction_elementary_transformation {
+ struct ga_instruction_elementary_trans_hess
+ : public ga_instruction_hess, ga_instruction_elementary_trans {
// Z(ndof,target_dim,N,N), coeff_in(Qmult,ndof) --> t(target_dim*Qmult,N,N)
virtual int exec() {
GA_DEBUG_INFO("Instruction: Hessian with elementary transformation");
- do_transformation();
+ size_type ndof = Z.sizes()[0];
+ size_type Qmult = qdim / Z.sizes()[1];
+ do_transformation(ndof*Qmult, t.sizes()[0]);
return ga_instruction_hess::exec();
}
- ga_instruction_elementary_transformation_hess
+ ga_instruction_elementary_trans_hess
(base_tensor &tt, const base_tensor &Z_, const base_vector &co, size_type
q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_hess(tt, Z_, coeff_out, q),
- ga_instruction_elementary_transformation(co, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans(co, e, mf1_, mf2_, ctx_, M_, icv_) {}
};
- struct ga_instruction_elementary_transformation_diverg
- : public ga_instruction_diverg, ga_instruction_elementary_transformation {
+ struct ga_instruction_elementary_trans_diverg
+ : public ga_instruction_diverg, ga_instruction_elementary_trans {
// Z(ndof,target_dim,N), coeff_in(Qmult,ndof) --> t(1)
virtual int exec() {
GA_DEBUG_INFO("Instruction: divergence with elementary transformation");
- do_transformation();
+ size_type ndof = Z.sizes()[0];
+ size_type Qmult = qdim / Z.sizes()[1];
+ do_transformation(ndof*Qmult, t.sizes()[0]);
return ga_instruction_diverg::exec();
}
- ga_instruction_elementary_transformation_diverg
+ ga_instruction_elementary_trans_diverg
(base_tensor &tt, const base_tensor &Z_, const base_vector &co, size_type
q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_diverg(tt, Z_, coeff_out, q),
- ga_instruction_elementary_transformation(co, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans(co, e, mf1_, mf2_, ctx_, M_, icv_) {}
};
struct ga_instruction_update_group_info : public ga_instruction {
@@ -1528,131 +1527,126 @@ namespace getfem {
};
- struct ga_instruction_elementary_transformation_base {
+ struct ga_instruction_elementary_trans_base {
base_tensor t_in;
base_tensor &t_out;
pelementary_transformation elemtrans;
- const mesh_fem &mf;
+ const mesh_fem &mf1, &mf2;
const fem_interpolation_context &ctx;
base_matrix &M;
- const mesh_fem **mf_M;
size_type &icv;
- void do_transformation(size_type n) {
- if (M.size() == 0 || icv != ctx.convex_num() || &mf != *mf_M) {
- M.base_resize(n, n);
- *mf_M = &mf; icv = ctx.convex_num();
- elemtrans->give_transformation(mf, icv, M);
+ void do_transformation(size_type n, size_type m) {
+ if (icv != ctx.convex_num() || M.size() == 0) {
+ M.base_resize(m, n);
+ icv = ctx.convex_num();
+ elemtrans->give_transformation(mf1, mf2, icv, M);
}
t_out.mat_reduction(t_in, M, 0);
}
- ga_instruction_elementary_transformation_base
- (base_tensor &t_, pelementary_transformation e, const mesh_fem &mf_,
- const fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
- : t_out(t_), elemtrans(e), mf(mf_), ctx(ctx_),
- M(M_), mf_M(mf_M_), icv(icv_) {}
+ ga_instruction_elementary_trans_base
+ (base_tensor &t_, pelementary_transformation e, const mesh_fem &mf1_,
+ const mesh_fem &mf2_,
+ const fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
+ : t_out(t_), elemtrans(e), mf1(mf1_), mf2(mf2_), ctx(ctx_),
+ M(M_), icv(icv_) {}
};
- struct ga_instruction_elementary_transformation_val_base
+ struct ga_instruction_elementary_trans_val_base
: public ga_instruction_copy_val_base,
- ga_instruction_elementary_transformation_base {
+ ga_instruction_elementary_trans_base {
// Z(ndof,target_dim) --> t_in --> t_out(Qmult*ndof,Qmult*target_dim)
virtual int exec() {
GA_DEBUG_INFO("Instruction: value of test functions with elementary "
"transformation");
size_type ndof = Z.sizes()[0];
size_type Qmult = qdim / Z.sizes()[1];
- t_in.adjust_sizes(t_out.sizes());
+ t_in.adjust_sizes(Qmult*ndof, Qmult*Z.sizes()[1]);
ga_instruction_copy_val_base::exec();
- do_transformation(ndof*Qmult);
+ do_transformation(ndof*Qmult, t_out.sizes()[0]);
return 0;
}
- ga_instruction_elementary_transformation_val_base
+ ga_instruction_elementary_trans_val_base
(base_tensor &t_, const base_tensor &Z_, size_type q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_copy_val_base(t_in, Z_, q),
- ga_instruction_elementary_transformation_base(t_, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans_base(t_, e, mf1_, mf2_, ctx_,
+ M_, icv_) {}
};
- struct ga_instruction_elementary_transformation_grad_base
+ struct ga_instruction_elementary_trans_grad_base
: public ga_instruction_copy_grad_base,
- ga_instruction_elementary_transformation_base {
+ ga_instruction_elementary_trans_base {
// Z(ndof,target_dim,N) --> t_in --> t_out(Qmult*ndof,Qmult*target_dim,N)
virtual int exec() {
GA_DEBUG_INFO("Instruction: gradient of test functions with elementary "
"transformation");
size_type ndof = Z.sizes()[0];
size_type Qmult = qdim / Z.sizes()[1];
- t_in.adjust_sizes(t_out.sizes());
+ t_in.adjust_sizes(Qmult*ndof, Qmult*Z.sizes()[1], Z.sizes()[2]);
ga_instruction_copy_grad_base::exec();
- do_transformation(ndof*Qmult);
+ do_transformation(ndof*Qmult, t_out.sizes()[0]);
return 0;
}
- ga_instruction_elementary_transformation_grad_base
+ ga_instruction_elementary_trans_grad_base
(base_tensor &t_, const base_tensor &Z_, size_type q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_copy_grad_base(t_in, Z_, q),
- ga_instruction_elementary_transformation_base(t_, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans_base(t_, e, mf1_, mf2_, ctx_,
+ M_, icv_) {}
};
- struct ga_instruction_elementary_transformation_hess_base
+ struct ga_instruction_elementary_trans_hess_base
: public ga_instruction_copy_hess_base,
- ga_instruction_elementary_transformation_base {
+ ga_instruction_elementary_trans_base {
// Z(ndof,target_dim,N*N) --> t_out(Qmult*ndof,Qmult*target_dim,N,N)
virtual int exec() {
GA_DEBUG_INFO("Instruction: Hessian of test functions with elementary "
"transformation");
size_type ndof = Z.sizes()[0];
size_type Qmult = qdim / Z.sizes()[1];
- t_in.adjust_sizes(t_out.sizes());
+ t_in.adjust_sizes(Qmult*ndof, Qmult*Z.sizes()[1], Z.sizes()[2]);
ga_instruction_copy_hess_base::exec();
- do_transformation(ndof*Qmult);
+ do_transformation(ndof*Qmult, t_out.sizes()[0]);
return 0;
}
- ga_instruction_elementary_transformation_hess_base
+ ga_instruction_elementary_trans_hess_base
(base_tensor &t_, const base_tensor &Z_, size_type q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_copy_hess_base(t_in, Z_, q),
- ga_instruction_elementary_transformation_base(t_, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans_base(t_, e, mf1_, mf2_, ctx_,
+ M_, icv_) {}
};
- struct ga_instruction_elementary_transformation_diverg_base
+ struct ga_instruction_elementary_trans_diverg_base
: public ga_instruction_copy_diverg_base,
- ga_instruction_elementary_transformation_base {
+ ga_instruction_elementary_trans_base {
// Z(ndof,target_dim,N) --> t_out(Qmult*ndof)
virtual int exec() {
GA_DEBUG_INFO("Instruction: divergence of test functions with elementary
"
"transformation");
size_type ndof = Z.sizes()[0];
size_type Qmult = qdim / Z.sizes()[1];
- t_in.adjust_sizes(t_out.sizes());
+ t_in.adjust_sizes(Qmult*ndof);
ga_instruction_copy_diverg_base::exec();
- do_transformation(ndof*Qmult);
+ do_transformation(ndof*Qmult, t_out.sizes()[0]);
return 0;
}
- ga_instruction_elementary_transformation_diverg_base
+ ga_instruction_elementary_trans_diverg_base
(base_tensor &t_, const base_tensor &Z_, size_type q,
- pelementary_transformation e, const mesh_fem &mf_,
- fem_interpolation_context &ctx_, base_matrix &M_,
- const mesh_fem **mf_M_, size_type &icv_)
+ pelementary_transformation e, const mesh_fem &mf1_, const mesh_fem &mf2_,
+ fem_interpolation_context &ctx_, base_matrix &M_, size_type &icv_)
: ga_instruction_copy_diverg_base(t_in, Z_, q),
- ga_instruction_elementary_transformation_base(t_, e, mf_, ctx_, M_,
- mf_M_, icv_) {}
+ ga_instruction_elementary_trans_base(t_, e, mf1_, mf2_, ctx_,
+ M_, icv_) {}
};
@@ -5063,287 +5057,310 @@ namespace getfem {
case GA_NODE_XFEM_PLUS_HESS: case GA_NODE_XFEM_PLUS_DIVERG:
case GA_NODE_XFEM_MINUS_VAL: case GA_NODE_XFEM_MINUS_GRAD:
case GA_NODE_XFEM_MINUS_HESS: case GA_NODE_XFEM_MINUS_DIVERG:
- if (function_case) {
- GMM_ASSERT1(pnode->node_type != GA_NODE_ELEMENTARY_VAL &&
- pnode->node_type != GA_NODE_ELEMENTARY_GRAD &&
- pnode->node_type != GA_NODE_ELEMENTARY_HESS &&
- pnode->node_type != GA_NODE_ELEMENTARY_DIVERG,
- "No elementary transformation is allowed in functions");
- GMM_ASSERT1(pnode->node_type != GA_NODE_XFEM_PLUS_VAL &&
- pnode->node_type != GA_NODE_XFEM_PLUS_GRAD &&
- pnode->node_type != GA_NODE_XFEM_PLUS_HESS &&
- pnode->node_type != GA_NODE_XFEM_PLUS_DIVERG,
- "Xfem_plus not allowed in functions");
- GMM_ASSERT1(pnode->node_type != GA_NODE_XFEM_MINUS_VAL &&
- pnode->node_type != GA_NODE_XFEM_MINUS_GRAD &&
- pnode->node_type != GA_NODE_XFEM_MINUS_HESS &&
- pnode->node_type != GA_NODE_XFEM_MINUS_DIVERG,
- "Xfem_plus not allowed in functions");
- const mesh_fem *mf = workspace.associated_mf(pnode->name);
- const im_data *imd = workspace.associated_im_data(pnode->name);
- GMM_ASSERT1(!mf,"No fem expression is allowed in function expression");
- GMM_ASSERT1(!imd, "No integration method data is allowed in "
- "function expression");
- if (gmm::vect_size(workspace.value(pnode->name)) == 1)
- pgai = std::make_shared<ga_instruction_copy_scalar>
- (pnode->tensor()[0], (workspace.value(pnode->name))[0]);
- else
- pgai = std::make_shared<ga_instruction_copy_vect>
- (pnode->tensor().as_vector(), workspace.value(pnode->name));
- rmi.instructions.push_back(std::move(pgai));
- } else {
- const mesh_fem *mf = workspace.associated_mf(pnode->name);
- const im_data *imd = workspace.associated_im_data(pnode->name);
-
- if (imd) {
- pgai = std::make_shared<ga_instruction_extract_local_im_data>
- (pnode->tensor(), *imd, workspace.value(pnode->name),
- gis.pai, gis.ctx, workspace.qdim(pnode->name));
+ {
+ bool is_elementary = (pnode->node_type == GA_NODE_ELEMENTARY_VAL ||
+ pnode->node_type == GA_NODE_ELEMENTARY_GRAD ||
+ pnode->node_type == GA_NODE_ELEMENTARY_HESS ||
+ pnode->node_type == GA_NODE_ELEMENTARY_DIVERG);
+ if (function_case) {
+ GMM_ASSERT1(!is_elementary,
+ "No elementary transformation is allowed in functions");
+ GMM_ASSERT1(pnode->node_type != GA_NODE_XFEM_PLUS_VAL &&
+ pnode->node_type != GA_NODE_XFEM_PLUS_GRAD &&
+ pnode->node_type != GA_NODE_XFEM_PLUS_HESS &&
+ pnode->node_type != GA_NODE_XFEM_PLUS_DIVERG,
+ "Xfem_plus not allowed in functions");
+ GMM_ASSERT1(pnode->node_type != GA_NODE_XFEM_MINUS_VAL &&
+ pnode->node_type != GA_NODE_XFEM_MINUS_GRAD &&
+ pnode->node_type != GA_NODE_XFEM_MINUS_HESS &&
+ pnode->node_type != GA_NODE_XFEM_MINUS_DIVERG,
+ "Xfem_plus not allowed in functions");
+ const mesh_fem *mf = workspace.associated_mf(pnode->name);
+ const im_data *imd = workspace.associated_im_data(pnode->name);
+ GMM_ASSERT1(!mf, "No fem expression is allowed in "
+ "function expression");
+ GMM_ASSERT1(!imd, "No integration method data is allowed in "
+ "function expression");
+ if (gmm::vect_size(workspace.value(pnode->name)) == 1)
+ pgai = std::make_shared<ga_instruction_copy_scalar>
+ (pnode->tensor()[0], (workspace.value(pnode->name))[0]);
+ else
+ pgai = std::make_shared<ga_instruction_copy_vect>
+ (pnode->tensor().as_vector(), workspace.value(pnode->name));
rmi.instructions.push_back(std::move(pgai));
} else {
- GMM_ASSERT1(mf, "Internal error");
-
- GMM_ASSERT1(&(mf->linked_mesh()) == &(m),
- "The finite element of variable " << pnode->name <<
- " has to be defined on the same mesh than the "
- "integration method or interpolation used");
-
- // An instruction for extracting local dofs of the variable.
- if (rmi.local_dofs.count(pnode->name) == 0) {
- rmi.local_dofs[pnode->name] = base_vector(1);
- extend_variable_in_gis(workspace, pnode->name, gis);
- // cout << "local dof of " << pnode->name << endl;
- size_type qmult2 = mf->get_qdim();
- if (qmult2 > 1 && !(mf->is_uniformly_vectorized()))
- qmult2 = size_type(-1);
- pgai = std::make_shared<ga_instruction_slice_local_dofs>
- (*mf, *(gis.extended_vars[pnode->name]), gis.ctx,
- rmi.local_dofs[pnode->name],
- workspace.qdim(pnode->name) / mf->get_qdim(), qmult2);
- rmi.elt_instructions.push_back(std::move(pgai));
+ const mesh_fem *mf = workspace.associated_mf(pnode->name), *mfo=mf;
+ const im_data *imd = workspace.associated_im_data(pnode->name);
+
+ if (is_elementary) {
+ mf = workspace.associated_mf(pnode->elementary_target);
+ GMM_ASSERT1(mf && mfo,
+ "Wrong context for elementary transformation");
+ GMM_ASSERT1(&(mfo->linked_mesh()) == &(m),
+ "The finite element of variable " << pnode->name
+ << " has to be defined on the same mesh than the "
+ << "integration method or interpolation used");
}
-
- // An instruction for pfp update
- if (mf->is_uniform()) {
- if (rmi.pfps.count(mf) == 0) {
+
+ if (imd) {
+ pgai = std::make_shared<ga_instruction_extract_local_im_data>
+ (pnode->tensor(), *imd, workspace.value(pnode->name),
+ gis.pai, gis.ctx, workspace.qdim(pnode->name));
+ rmi.instructions.push_back(std::move(pgai));
+ } else {
+ GMM_ASSERT1(mf, "Internal error");
+
+ GMM_ASSERT1(&(mf->linked_mesh()) == &(m),
+ "The finite element of variable " <<
+ (is_elementary ? pnode->elementary_target :
pnode->name)
+ << " has to be defined on the same mesh than the "
+ << "integration method or interpolation used");
+
+ // An instruction for extracting local dofs of the variable.
+ if (rmi.local_dofs.count(pnode->name) == 0) {
+ rmi.local_dofs[pnode->name] = base_vector(1);
+ extend_variable_in_gis(workspace, pnode->name, gis);
+ // cout << "local dof of " << pnode->name << endl;
+ size_type qmult2 = mfo->get_qdim();
+ if (qmult2 > 1 && !(mfo->is_uniformly_vectorized()))
+ qmult2 = size_type(-1);
+ pgai = std::make_shared<ga_instruction_slice_local_dofs>
+ (*mfo, *(gis.extended_vars[pnode->name]), gis.ctx,
+ rmi.local_dofs[pnode->name],
+ workspace.qdim(pnode->name) / mfo->get_qdim(), qmult2);
+ rmi.elt_instructions.push_back(std::move(pgai));
+ }
+
+ // An instruction for pfp update
+ if (mf->is_uniform()) {
+ if (rmi.pfps.count(mf) == 0) {
+ rmi.pfps[mf] = 0;
+ pgai = std::make_shared<ga_instruction_update_pfp>
+ (*mf, rmi.pfps[mf], gis.ctx, gis.fp_pool);
+ rmi.begin_instructions.push_back(std::move(pgai));
+ }
+ } else if (rmi.pfps.count(mf) == 0 ||
+ !if_hierarchy.is_compatible(rmi.pfp_hierarchy[mf])) {
+ rmi.pfp_hierarchy[mf].push_back(if_hierarchy);
rmi.pfps[mf] = 0;
pgai = std::make_shared<ga_instruction_update_pfp>
(*mf, rmi.pfps[mf], gis.ctx, gis.fp_pool);
- rmi.begin_instructions.push_back(std::move(pgai));
- }
- } else if (rmi.pfps.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.pfp_hierarchy[mf])) {
- rmi.pfp_hierarchy[mf].push_back(if_hierarchy);
- rmi.pfps[mf] = 0;
- pgai = std::make_shared<ga_instruction_update_pfp>
- (*mf, rmi.pfps[mf], gis.ctx, gis.fp_pool);
- rmi.instructions.push_back(std::move(pgai));
- }
-
- // An instruction for the base value
- pgai = pga_instruction();
- switch (pnode->node_type) {
- case GA_NODE_VAL: case GA_NODE_ELEMENTARY_VAL:
- if (rmi.base.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.base_hierarchy[mf])) {
- rmi.base_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_val_base>
- (rmi.base[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
- case GA_NODE_XFEM_PLUS_VAL:
- if (rmi.xfem_plus_base.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.xfem_plus_base_hierarchy[mf]))
- {
- rmi.xfem_plus_base_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_xfem_plus_val_base>
- (rmi.xfem_plus_base[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
- case GA_NODE_XFEM_MINUS_VAL:
- if (rmi.xfem_minus_base.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.xfem_minus_base_hierarchy[mf]))
- {
- rmi.xfem_minus_base_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_xfem_minus_val_base>
- (rmi.xfem_minus_base[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
- case GA_NODE_GRAD: case GA_NODE_DIVERG:
- case GA_NODE_ELEMENTARY_GRAD: case GA_NODE_ELEMENTARY_DIVERG:
- if (rmi.grad.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.grad_hierarchy[mf])) {
- rmi.grad_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_grad_base>
- (rmi.grad[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
- case GA_NODE_XFEM_PLUS_GRAD: case GA_NODE_XFEM_PLUS_DIVERG:
- if (rmi.xfem_plus_grad.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.xfem_plus_grad_hierarchy[mf]))
- {
- rmi.xfem_plus_grad_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_xfem_plus_grad_base>
- (rmi.xfem_plus_grad[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
- case GA_NODE_XFEM_MINUS_GRAD: case GA_NODE_XFEM_MINUS_DIVERG:
- if (rmi.xfem_minus_grad.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.xfem_minus_grad_hierarchy[mf]))
- {
- rmi.xfem_minus_grad_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_xfem_minus_grad_base>
- (rmi.xfem_minus_grad[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ rmi.instructions.push_back(std::move(pgai));
}
- break;
- case GA_NODE_HESS: case GA_NODE_ELEMENTARY_HESS:
- if (rmi.hess.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.hess_hierarchy[mf])) {
- rmi.hess_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_hess_base>
- (rmi.hess[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
- case GA_NODE_XFEM_PLUS_HESS:
- if (rmi.xfem_plus_hess.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.xfem_plus_hess_hierarchy[mf]))
- {
- rmi.xfem_plus_hess_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_xfem_plus_hess_base>
- (rmi.xfem_plus_hess[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
- case GA_NODE_XFEM_MINUS_HESS:
- if (rmi.xfem_minus_hess.count(mf) == 0 ||
- !if_hierarchy.is_compatible(rmi.xfem_minus_hess_hierarchy[mf]))
- {
- rmi.xfem_minus_hess_hierarchy[mf].push_back(if_hierarchy);
- pgai = std::make_shared<ga_instruction_xfem_minus_hess_base>
- (rmi.xfem_minus_hess[mf], gis.ctx, *mf, rmi.pfps[mf]);
- }
- break;
-
- default : GMM_ASSERT1(false, "Internal error");
- }
- if (pgai) rmi.instructions.push_back(std::move(pgai));
-
- // The eval instruction
- switch (pnode->node_type) {
- case GA_NODE_VAL: // --> t(target_dim*Qmult)
- pgai = std::make_shared<ga_instruction_val>
- (pnode->tensor(), rmi.base[mf], rmi.local_dofs[pnode->name],
- workspace.qdim(pnode->name));
- break;
- case GA_NODE_GRAD: // --> t(target_dim*Qmult,N)
- pgai = std::make_shared<ga_instruction_grad>
- (pnode->tensor(), rmi.grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_HESS: // --> t(target_dim*Qmult,N,N)
- pgai = std::make_shared<ga_instruction_hess>
- (pnode->tensor(), rmi.hess[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_DIVERG: // --> t(1)
- pgai = std::make_shared<ga_instruction_diverg>
- (pnode->tensor(), rmi.grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_PLUS_VAL: // --> t(target_dim*Qmult)
- pgai = std::make_shared<ga_instruction_val>
- (pnode->tensor(), rmi.xfem_plus_base[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_PLUS_GRAD: // --> t(target_dim*Qmult,N)
- pgai = std::make_shared<ga_instruction_grad>
- (pnode->tensor(), rmi.xfem_plus_grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_PLUS_HESS: // --> t(target_dim*Qmult,N,N)
- pgai = std::make_shared<ga_instruction_hess>
- (pnode->tensor(), rmi.xfem_plus_hess[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_PLUS_DIVERG: // --> t(1)
- pgai = std::make_shared<ga_instruction_diverg>
- (pnode->tensor(), rmi.xfem_plus_grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_MINUS_VAL: // --> t(target_dim*Qmult)
- pgai = std::make_shared<ga_instruction_val>
- (pnode->tensor(), rmi.xfem_minus_base[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_MINUS_GRAD: // --> t(target_dim*Qmult,N)
- pgai = std::make_shared<ga_instruction_grad>
- (pnode->tensor(), rmi.xfem_minus_grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_MINUS_HESS: // --> t(target_dim*Qmult,N,N)
- pgai = std::make_shared<ga_instruction_hess>
- (pnode->tensor(), rmi.xfem_minus_hess[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_XFEM_MINUS_DIVERG: // --> t(1)
- pgai = std::make_shared<ga_instruction_diverg>
- (pnode->tensor(), rmi.xfem_minus_grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
- break;
- case GA_NODE_ELEMENTARY_VAL:
- { // --> t(target_dim*Qmult)
- ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
- pgai =
- std::make_shared<ga_instruction_elementary_transformation_val>
- (pnode->tensor(), rmi.base[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name),
- workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
+
+ // An instruction for the base value
+ pgai = pga_instruction();
+ switch (pnode->node_type) {
+ case GA_NODE_VAL: case GA_NODE_ELEMENTARY_VAL:
+ if (rmi.base.count(mf) == 0 ||
+ !if_hierarchy.is_compatible(rmi.base_hierarchy[mf])) {
+ rmi.base_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_val_base>
+ (rmi.base[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_XFEM_PLUS_VAL:
+ if (rmi.xfem_plus_base.count(mf) == 0 ||
+
!if_hierarchy.is_compatible(rmi.xfem_plus_base_hierarchy[mf]))
+ {
+ rmi.xfem_plus_base_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_xfem_plus_val_base>
+ (rmi.xfem_plus_base[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_XFEM_MINUS_VAL:
+ if (rmi.xfem_minus_base.count(mf) == 0 ||
+
!if_hierarchy.is_compatible(rmi.xfem_minus_base_hierarchy[mf]))
+ {
+ rmi.xfem_minus_base_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_xfem_minus_val_base>
+ (rmi.xfem_minus_base[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_GRAD: case GA_NODE_DIVERG:
+ case GA_NODE_ELEMENTARY_GRAD: case GA_NODE_ELEMENTARY_DIVERG:
+ if (rmi.grad.count(mf) == 0 ||
+ !if_hierarchy.is_compatible(rmi.grad_hierarchy[mf])) {
+ rmi.grad_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_grad_base>
+ (rmi.grad[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_XFEM_PLUS_GRAD: case GA_NODE_XFEM_PLUS_DIVERG:
+ if (rmi.xfem_plus_grad.count(mf) == 0 ||
+
!if_hierarchy.is_compatible(rmi.xfem_plus_grad_hierarchy[mf]))
+ {
+ rmi.xfem_plus_grad_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_xfem_plus_grad_base>
+ (rmi.xfem_plus_grad[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_XFEM_MINUS_GRAD: case GA_NODE_XFEM_MINUS_DIVERG:
+ if (rmi.xfem_minus_grad.count(mf) == 0 ||
+
!if_hierarchy.is_compatible(rmi.xfem_minus_grad_hierarchy[mf]))
+ {
+ rmi.xfem_minus_grad_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_xfem_minus_grad_base>
+ (rmi.xfem_minus_grad[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_HESS: case GA_NODE_ELEMENTARY_HESS:
+ if (rmi.hess.count(mf) == 0 ||
+ !if_hierarchy.is_compatible(rmi.hess_hierarchy[mf])) {
+ rmi.hess_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_hess_base>
+ (rmi.hess[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_XFEM_PLUS_HESS:
+ if (rmi.xfem_plus_hess.count(mf) == 0 ||
+
!if_hierarchy.is_compatible(rmi.xfem_plus_hess_hierarchy[mf]))
+ {
+ rmi.xfem_plus_hess_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_xfem_plus_hess_base>
+ (rmi.xfem_plus_hess[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+ case GA_NODE_XFEM_MINUS_HESS:
+ if (rmi.xfem_minus_hess.count(mf) == 0 ||
+
!if_hierarchy.is_compatible(rmi.xfem_minus_hess_hierarchy[mf]))
+ {
+ rmi.xfem_minus_hess_hierarchy[mf].push_back(if_hierarchy);
+ pgai = std::make_shared<ga_instruction_xfem_minus_hess_base>
+ (rmi.xfem_minus_hess[mf], gis.ctx, *mf, rmi.pfps[mf]);
+ }
+ break;
+
+ default : GMM_ASSERT1(false, "Internal error");
}
- break;
- case GA_NODE_ELEMENTARY_GRAD:
- { // --> t(target_dim*Qmult,N)
- ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
- pgai =
- std::make_shared<ga_instruction_elementary_transformation_grad>
+ if (pgai) rmi.instructions.push_back(std::move(pgai));
+
+ // The eval instruction
+ switch (pnode->node_type) {
+ case GA_NODE_VAL: // --> t(target_dim*Qmult)
+ pgai = std::make_shared<ga_instruction_val>
+ (pnode->tensor(), rmi.base[mf], rmi.local_dofs[pnode->name],
+ workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_GRAD: // --> t(target_dim*Qmult,N)
+ pgai = std::make_shared<ga_instruction_grad>
(pnode->tensor(), rmi.grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name),
- workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
- }
- break;
- case GA_NODE_ELEMENTARY_HESS:
- { // --> t(target_dim*Qmult,N,N)
- ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
- pgai =
- std::make_shared<ga_instruction_elementary_transformation_hess>
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_HESS: // --> t(target_dim*Qmult,N,N)
+ pgai = std::make_shared<ga_instruction_hess>
(pnode->tensor(), rmi.hess[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name),
- workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
- }
- break;
- case GA_NODE_ELEMENTARY_DIVERG:
- { // --> t(1)
- ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
- pgai =
-
std::make_shared<ga_instruction_elementary_transformation_diverg>
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_DIVERG: // --> t(1)
+ pgai = std::make_shared<ga_instruction_diverg>
(pnode->tensor(), rmi.grad[mf],
- rmi.local_dofs[pnode->name], workspace.qdim(pnode->name),
- workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_PLUS_VAL: // --> t(target_dim*Qmult)
+ pgai = std::make_shared<ga_instruction_val>
+ (pnode->tensor(), rmi.xfem_plus_base[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_PLUS_GRAD: // --> t(target_dim*Qmult,N)
+ pgai = std::make_shared<ga_instruction_grad>
+ (pnode->tensor(), rmi.xfem_plus_grad[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_PLUS_HESS: // --> t(target_dim*Qmult,N,N)
+ pgai = std::make_shared<ga_instruction_hess>
+ (pnode->tensor(), rmi.xfem_plus_hess[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_PLUS_DIVERG: // --> t(1)
+ pgai = std::make_shared<ga_instruction_diverg>
+ (pnode->tensor(), rmi.xfem_plus_grad[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_MINUS_VAL: // --> t(target_dim*Qmult)
+ pgai = std::make_shared<ga_instruction_val>
+ (pnode->tensor(), rmi.xfem_minus_base[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_MINUS_GRAD: // --> t(target_dim*Qmult,N)
+ pgai = std::make_shared<ga_instruction_grad>
+ (pnode->tensor(), rmi.xfem_minus_grad[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_MINUS_HESS: // --> t(target_dim*Qmult,N,N)
+ pgai = std::make_shared<ga_instruction_hess>
+ (pnode->tensor(), rmi.xfem_minus_hess[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_XFEM_MINUS_DIVERG: // --> t(1)
+ pgai = std::make_shared<ga_instruction_diverg>
+ (pnode->tensor(), rmi.xfem_minus_grad[mf],
+ rmi.local_dofs[pnode->name], workspace.qdim(pnode->name));
+ break;
+ case GA_NODE_ELEMENTARY_VAL:
+ { // --> t(target_dim*Qmult)
+ ga_instruction_set::elementary_trans_info &eti
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
+ pgai =
+ std::make_shared<ga_instruction_elementary_trans_val>
+ (pnode->tensor(), rmi.base[mf],
+ rmi.local_dofs[pnode->name],
+ workspace.qdim(pnode->elementary_target),
+ workspace.elementary_transformation(pnode->elementary_name),
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
+ }
+ break;
+ case GA_NODE_ELEMENTARY_GRAD:
+ { // --> t(target_dim*Qmult,N)
+ ga_instruction_set::elementary_trans_info &eti
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
+ pgai =
+ std::make_shared<ga_instruction_elementary_trans_grad>
+ (pnode->tensor(), rmi.grad[mf],
+ rmi.local_dofs[pnode->name],
+ workspace.qdim(pnode->elementary_target),
+ workspace.elementary_transformation(pnode->elementary_name),
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
+ }
+ break;
+ case GA_NODE_ELEMENTARY_HESS:
+ { // --> t(target_dim*Qmult,N,N)
+ ga_instruction_set::elementary_trans_info &eti
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
+ pgai =
+ std::make_shared<ga_instruction_elementary_trans_hess>
+ (pnode->tensor(), rmi.hess[mf],
+ rmi.local_dofs[pnode->name],
+ workspace.qdim(pnode->elementary_target),
+ workspace.elementary_transformation(pnode->elementary_name),
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
+ }
+ break;
+ case GA_NODE_ELEMENTARY_DIVERG:
+ { // --> t(1)
+ ga_instruction_set::elementary_trans_info &eti
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
+ pgai =
+ std::make_shared<ga_instruction_elementary_trans_diverg>
+ (pnode->tensor(), rmi.grad[mf],
+ rmi.local_dofs[pnode->name],
+ workspace.qdim(pnode->elementary_target),
+ workspace.elementary_transformation(pnode->elementary_name),
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
+ }
+ break;
+ default: break;
}
- break;
- default: break;
+ rmi.instructions.push_back(std::move(pgai));
}
- rmi.instructions.push_back(std::move(pgai));
}
}
break;
-
+
case GA_NODE_SECONDARY_DOMAIN_VAL: case GA_NODE_SECONDARY_DOMAIN_GRAD:
case GA_NODE_SECONDARY_DOMAIN_HESS: case GA_NODE_SECONDARY_DOMAIN_DIVERG:
{
@@ -5527,12 +5544,27 @@ namespace getfem {
// GMM_ASSERT1(!function_case,
// "Test functions not allowed in functions");
{
- const mesh_fem *mf = workspace.associated_mf(pnode->name);
+ bool is_elementary = (pnode->node_type==GA_NODE_ELEMENTARY_VAL_TEST ||
+ pnode->node_type==GA_NODE_ELEMENTARY_GRAD_TEST ||
+ pnode->node_type==GA_NODE_ELEMENTARY_HESS_TEST ||
+
pnode->node_type==GA_NODE_ELEMENTARY_DIVERG_TEST);
+ const mesh_fem *mf = workspace.associated_mf(pnode->name), *mfo=mf;
+ if (is_elementary) {
+ mf = workspace.associated_mf(pnode->elementary_target);
+ GMM_ASSERT1(mf && mfo,
+ "Wrong context for elementary transformation");
+ GMM_ASSERT1(&(mfo->linked_mesh()) == &(m),
+ "The finite element of variable " << pnode->name
+ << " has to be defined on the same mesh than the "
+ << "integration method or interpolation used");
+ }
+
if (mf) {
GMM_ASSERT1(&(mf->linked_mesh()) == &(m),
"The finite element of variable " << pnode->name <<
- " and the applied integration method have to be"
- " defined on the same mesh");
+ (is_elementary ? pnode->elementary_target : pnode->name)
+ << " and the applied integration method have to be"
+ << " defined on the same mesh");
// An instruction for pfp update
if (is_uniform) {
@@ -5733,45 +5765,49 @@ namespace getfem {
case GA_NODE_ELEMENTARY_VAL_TEST:
{ // --> t(Qmult*ndof,Qmult*target_dim)
ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
pgai =
-
std::make_shared<ga_instruction_elementary_transformation_val_base>
+ std::make_shared<ga_instruction_elementary_trans_val_base>
(pnode->tensor(), rmi.base[mf], mf->get_qdim(),
workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
}
break;
case GA_NODE_ELEMENTARY_GRAD_TEST:
{ // --> t(Qmult*ndof,Qmult*target_dim,N)
ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
pgai =
-
std::make_shared<ga_instruction_elementary_transformation_grad_base>
+ std::make_shared<ga_instruction_elementary_trans_grad_base>
(pnode->tensor(), rmi.grad[mf], mf->get_qdim(),
workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
}
break;
case GA_NODE_ELEMENTARY_HESS_TEST:
{ // --> t(Qmult*ndof,Qmult*target_dim,N,N)
ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
pgai =
-
std::make_shared<ga_instruction_elementary_transformation_hess_base>
+ std::make_shared<ga_instruction_elementary_trans_hess_base>
(pnode->tensor(), rmi.hess[mf], mf->get_qdim(),
workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
}
break;
case GA_NODE_ELEMENTARY_DIVERG_TEST:
{ // --> t(Qmult*ndof)
ga_instruction_set::elementary_trans_info &eti
- = rmi.elementary_trans_infos[pnode->elementary_name];
+ = rmi.elementary_trans_infos
+ [std::make_tuple(pnode->elementary_name, mfo, mf)];
pgai =
-
std::make_shared<ga_instruction_elementary_transformation_diverg_base>
+ std::make_shared<ga_instruction_elementary_trans_diverg_base>
(pnode->tensor(), rmi.grad[mf], mf->get_qdim(),
workspace.elementary_transformation(pnode->elementary_name),
- *mf, gis.ctx, eti.M, &(eti.mf), eti.icv);
+ *mfo, *mf, gis.ctx, eti.M, eti.icv);
}
break;
default: break;
diff --git a/src/getfem_generic_assembly_semantic.cc
b/src/getfem_generic_assembly_semantic.cc
index cf27450..719d4da 100644
--- a/src/getfem_generic_assembly_semantic.cc
+++ b/src/getfem_generic_assembly_semantic.cc
@@ -1,6 +1,6 @@
/*===========================================================================
- Copyright (C) 2013-2018 Yves Renard
+ Copyright (C) 2013-2019 Yves Renard
This file is a part of GetFEM++
@@ -349,7 +349,8 @@ namespace getfem {
case GA_NODE_ELEMENTARY_VAL_TEST: case GA_NODE_ELEMENTARY_GRAD_TEST:
case GA_NODE_ELEMENTARY_HESS_TEST: case GA_NODE_ELEMENTARY_DIVERG_TEST:
c += 1.33*(1.22+ga_hash_code(pnode->name))
- + 2.63*ga_hash_code(pnode->elementary_name);
+ + 2.63*ga_hash_code(pnode->elementary_name)
+ + 3.47*ga_hash_code(pnode->elementary_target);
break;
case GA_NODE_XFEM_PLUS_VAL: case GA_NODE_XFEM_PLUS_GRAD:
case GA_NODE_XFEM_PLUS_HESS: case GA_NODE_XFEM_PLUS_DIVERG:
@@ -585,6 +586,13 @@ namespace getfem {
size_type test = ga_parse_prefix_test(name);
pnode->name = name;
+ if (ndt == 2) {
+ std::string target_name = pnode->elementary_target;
+ ga_parse_prefix_operator(target_name);
+ ga_parse_prefix_test(target_name);
+ pnode->elementary_target = target_name;
+ }
+
// Group must be tested and it should be a fem variable
if (!(workspace.variable_or_group_exists(name)))
ga_throw_error(pnode->expr, pnode->pos,
@@ -763,6 +771,15 @@ namespace getfem {
ga_throw_error(pnode->expr, pnode->pos,
"Unknown elementary transformation");
}
+ if (!(workspace.variable_or_group_exists(pnode->elementary_target)))
{
+ ga_throw_error(pnode->expr, pnode->pos, "Unknown data or variable "
+ << pnode->elementary_target);
+ }
+ const mesh_fem *mft = workspace.associated_mf(name);
+ if (!mft)
+ ga_throw_error(pnode->expr, pnode->pos,
+ "Thir argument of the elementary transformation "
+ "should be a finite element variables/data");
} else if (ndt == 3) {
if (!(workspace.secondary_domain_exists
(pnode->interpolate_name))) {
diff --git a/src/getfem_generic_assembly_tree.cc
b/src/getfem_generic_assembly_tree.cc
index 3617964..a9d1364 100644
--- a/src/getfem_generic_assembly_tree.cc
+++ b/src/getfem_generic_assembly_tree.cc
@@ -1,6 +1,6 @@
/*===========================================================================
- Copyright (C) 2013-2018 Yves Renard
+ Copyright (C) 2013-2019 Yves Renard
This file is a part of GetFEM++
@@ -1216,9 +1216,12 @@ namespace getfem {
if (is_interpolate)
str << "," << pnode->interpolate_name << ")";
- else if (is_elementary)
- str << "," << pnode->elementary_name << ")";
- else if (is_secondary)
+ else if (is_elementary) {
+ str << "," << pnode->elementary_name;
+ if (pnode->name.compare(pnode->elementary_target) != 0)
+ str << "," << pnode->elementary_target;
+ str << ")";
+ } else if (is_secondary)
str << ")";
else if (is_xfem_plus || is_xfem_minus)
str << ")";
@@ -1738,6 +1741,7 @@ namespace getfem {
"should be a variable or a test function.");
tree.current_node->name = std::string(&((*expr)[token_pos]),
token_length);
+ tree.current_node->elementary_target = tree.current_node->name;
t_type = ga_get_token(*expr, pos, token_pos, token_length);
if (t_type != GA_COMMA)
@@ -1751,6 +1755,19 @@ namespace getfem {
tree.current_node->elementary_name
= std::string(&((*expr)[token_pos]), token_length);
t_type = ga_get_token(*expr, pos, token_pos, token_length);
+
+ if (t_type == GA_COMMA) {
+ t_type = ga_get_token(*expr, pos, token_pos, token_length);
+ if (t_type != GA_NAME)
+ ga_throw_error(expr, pos,
+ "Third argument of Elementary_transformation "
+ "should be a variable or data name.");
+
+ tree.current_node->elementary_target =
+ std::string(&((*expr)[token_pos]), token_length);
+ t_type = ga_get_token(*expr, pos, token_pos, token_length);
+ }
+
if (t_type != GA_RPAR)
ga_throw_error(expr, pos-1, "Missing a parenthesis after "
"Elementary_transformation arguments.");
diff --git a/src/getfem_linearized_plates.cc b/src/getfem_linearized_plates.cc
index 3252757..a984f35 100644
--- a/src/getfem_linearized_plates.cc
+++ b/src/getfem_linearized_plates.cc
@@ -195,12 +195,15 @@ namespace getfem {
public:
- virtual void give_transformation(const mesh_fem &mf, size_type cv,
- base_matrix &M) const{
+ virtual void give_transformation(const mesh_fem &mf, const mesh_fem &mf2,
+ size_type cv, base_matrix &M) const{
THREAD_SAFE_STATIC base_matrix M_old;
THREAD_SAFE_STATIC pfem pf_old = nullptr;
-
+
+ GMM_ASSERT1(&mf == &mf2,
+ "This transformation works on identical fems only");
+
// Obtaining the fem descriptors
pfem pf1 = mf.fem_of_element(cv);
size_type N = 2;
@@ -296,11 +299,14 @@ namespace getfem {
public:
- virtual void give_transformation(const mesh_fem &mf, size_type cv,
- base_matrix &M) const{
+ virtual void give_transformation(const mesh_fem &mf, const mesh_fem &mf2,
+ size_type cv, base_matrix &M) const{
THREAD_SAFE_STATIC base_matrix M_old;
THREAD_SAFE_STATIC pfem pf_old = nullptr;
+
+ GMM_ASSERT1(&mf == &mf2,
+ "This transformation works on identical fems only");
// Obtaining the fem descriptors
pfem pf1 = mf.fem_of_element(cv);
- [Getfem-commits] [getfem-commits] master updated (860f866 -> e920df7), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject), Yves Renard, 2019/08/27
- [Getfem-commits] (no subject),
Yves Renard <=