v0.12.1
BasicFeTools.hpp
Go to the documentation of this file.
1 /* This file is part of MoFEM.
2  * MoFEM is free software: you can redistribute it and/or modify it under
3  * the terms of the GNU Lesser General Public License as published by the
4  * Free Software Foundation, either version 3 of the License, or (at your
5  * option) any later version.
6  *
7  * MoFEM is distributed in the hope that it will be useful, but WITHOUT
8  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
9  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
10  * License for more details.
11  *
12  * You should have received a copy of the GNU Lesser General Public
13  * License along with MoFEM. If not, see <http://www.gnu.org/licenses/>. */
14 // #include <BasicFiniteElements.hpp>
15 // namespace ContactPlasticity {
16 #pragma once
17 
18 struct CacheParams {
19 
20  double young_modulus;
22  double poisson;
23  double density;
24  // for plastic
25  double sigmaY;
26  double cn_pl;
27  // hardening
28  double H;
29  double C1_k;
30  double Q_inf;
31  double b_iso;
32  double visH;
34 
35  MatrixDouble mtD; // for stiffness matrix
36 
37  // for contact
38  double cn_cont;
40 
41  // for rotation
42  static vector<double> gravity;
43  static double spring_stiffness;
44 
45  static int closest_roller_id;
46  static double delta;
48  static Tensor2<double, 3, 3> Omega;
49  static Tensor4<double, 3, 3, 3, 3> Is; // symmetric
50 
51  static boost::ptr_vector<RigidBodyData> rollerDataVec;
53 
56  sigmaY(1e9), cn_pl(1), H(0), C1_k(0), Q_inf(0), b_iso(0), visH(0),
57  cn_cont(1. / young_modulus), /*spring_stiffness(0),
58  closest_roller_id(0),*/
59  rollers_stop_time(INFINITY), scale_constraint(1) {
60  mtD.resize(36, 1);
61  }
62 
64 
65  template <typename T1, typename T2>
67 
68  std::vector<double> gaps;
70  for (auto &roller : rollerDataVec) {
71  auto t_normal = roller.getNormal(t_coords, t_disp);
72  double dist = roller.getGap(t_coords);
73  gaps.push_back(dist);
74  }
75 
76  auto it = std::min_element(gaps.begin(), gaps.end());
77  auto nb = distance(gaps.begin(), it);
78 
79  return nb;
80  };
81 };
82 
91 
98 
99 static Number<0> N0;
100 static Number<1> N1;
101 static Number<2> N2;
102 
103 struct BlockParamData : public CacheParams {
104 
107 
108  // CHKERR PetscOptionsInsertString(NULL, default_options);
109  CHKERR PetscOptionsInsertString(
110  NULL,
111  "-mat_mumps_icntl_14 1200 -mat_mumps_icntl_24 1 -mat_mumps_icntl_13 1 ");
112 
113  char mumps_options[] = "-fieldsplit_0_mat_mumps_icntl_14 800 "
114  "-fieldsplit_0_mat_mumps_icntl_24 1 "
115  "-fieldsplit_0_mat_mumps_icntl_13 1 "
116  "-fieldsplit_1_mat_mumps_icntl_14 800 "
117  "-fieldsplit_1_mat_mumps_icntl_24 1 "
118  "-fieldsplit_1_mat_mumps_icntl_13 1";
119  CHKERR PetscOptionsInsertString(NULL, mumps_options);
120 
121  CHKERR PetscOptionsBegin(PETSC_COMM_WORLD, "", "", "none");
122 
123  int nb_ent_grav = 3;
124  CHKERR PetscOptionsGetRealArray(NULL, NULL, "-gravity", gravity.data(),
125  &nb_ent_grav, PETSC_NULL);
126  // for elastic
127  CHKERR PetscOptionsScalar("-young", "Young's modulus", "", young_modulus,
128  &young_modulus, PETSC_NULL);
129 
130  CHKERR PetscOptionsScalar("-poisson", "Poisson ratio", "", poisson,
131  &poisson, PETSC_NULL);
132  CHKERR PetscOptionsScalar("-density", "Density", "", density, &density,
133  PETSC_NULL);
134 
135  PetscBool with_contact = PETSC_FALSE;
136  CHKERR PetscOptionsBool("-with_contact", "run calculations with contact",
137  "", with_contact, &with_contact, PETSC_NULL);
138 
139  double rot_vec[3] = {0, 0, 1};
140  int nb_ent = 3;
141  PetscBool is_rotating;
142  // get rotation velocity vector
143  CHKERR PetscOptionsGetRealArray(NULL, NULL, "-rot_omega", rot_vec, &nb_ent,
144  &is_rotating);
145  if (nb_ent < 3 && is_rotating) {
146  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
147  "provide an appropriate number of entries for omega vector (3)");
148  }
149  {
150  Tensor1<double, 3> t1_omega(rot_vec[0], rot_vec[1], rot_vec[2]);
151  omega(i) = t1_omega(i);
152  Omega(i, k) = levi_civita<double>(i, j, k) * t1_omega(j);
153  }
154 
155  if (poisson >= 0.5)
156  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
157  "poisson ratio value not supported");
158 
159  CHKERR PetscOptionsScalar("-sigmaY", "SigmaY - limit stress", "", sigmaY,
160  &sigmaY, PETSC_NULL);
161  CHKERR PetscOptionsScalar("-H", "Hardening modulus", "", H, &H, PETSC_NULL);
162  CHKERR PetscOptionsScalar("-cn_pl",
163  "cn_pl parameter for active set constrains", "",
164  cn_pl, &cn_pl, PETSC_NULL);
165  CHKERR PetscOptionsScalar("-scale_constraint", "constraint scale", "",
166  scale_constraint, &scale_constraint, PETSC_NULL);
167  CHKERR PetscOptionsScalar("-C1_k", "1st Backstress", "", C1_k, &C1_k,
168  PETSC_NULL);
169  CHKERR PetscOptionsScalar("-Q_inf", "Saturation stress", "", Q_inf, &Q_inf,
170  PETSC_NULL);
171  CHKERR PetscOptionsScalar("-b_iso", "Isotropic exponent", "", b_iso, &b_iso,
172  PETSC_NULL);
173  CHKERR PetscOptionsScalar("-visH", "viscous hardening", "", visH, &visH,
174  PETSC_NULL);
175  // for contact
176  CHKERR PetscOptionsScalar("-spring",
177  "Springs stiffness on the boundary (contact)", "",
178  spring_stiffness, &spring_stiffness, PETSC_NULL);
179 
180  CHKERR PetscOptionsScalar("-cn_cont",
181  "cn_cont parameter for active set constrains", "",
182  cn_cont, &cn_cont, PETSC_NULL);
183  CHKERR PetscOptionsScalar(
184  "-rollers_stop", "Time at which rollers will stop moving", "",
185  rollers_stop_time, &rollers_stop_time, PETSC_NULL);
186 
187  constexpr int MAX_NB_OF_ROLLERS = 10;
188  // if (with_contact) {
189 
190  const char *body_list[] = {"plane", "sphere", "cylinder", "cone",
191  "torus", "roller", "stl", "nurbs"};
192  for (int ii = 0; ii < MAX_NB_OF_ROLLERS; ++ii) {
193  int nb_coord = 3;
194  int nb_disp = 3;
195  int nb_dirs = 3;
196  VectorDouble center_coords({0, 0, 0});
197  VectorDouble disp({0, 0, 0});
198  VectorDouble direction({0, 0, 1});
199  PetscBool flg = PETSC_FALSE;
200  PetscBool rflg;
201  string param = "-body" + to_string(ii);
202  // char mesh_file_name[255];
203  // CHKERR PetscOptionsString(param.c_str(), "file name", "", "mesh.vtk",
204  // mesh_file_name, 255, &flg);
205  int body_type_id = LASTBODY;
206  CHKERR PetscOptionsGetEList(PETSC_NULL, NULL, param.c_str(), body_list,
207  LASTBODY, &body_type_id, &flg);
208 
209  param = "-coords" + to_string(ii);
210  CHKERR PetscOptionsGetRealArray(NULL, NULL, param.c_str(),
211  &center_coords(0), &nb_coord, &rflg);
212  param = "-move" + to_string(ii);
213  CHKERR PetscOptionsGetRealArray(NULL, NULL, param.c_str(), &disp(0),
214  &nb_disp, &rflg);
215  if (flg) {
216 
217  switch (body_type_id) {
218  case PLANE:
219  rollerDataVec.push_back(new PlaneRigidBody(center_coords, disp, ii));
220  break;
221  case SPHERE:
222  rollerDataVec.push_back(new SphereRigidBody(center_coords, disp, ii));
223  break;
224  case CYLINDER:
225  rollerDataVec.push_back(
226  new CylinderRigidBody(center_coords, disp, ii));
227  break;
228  case CONE:
229  rollerDataVec.push_back(new ConeRigidBody(center_coords, disp, ii));
230  break;
231  case TORUS:
232  rollerDataVec.push_back(new TorusRigidBody(center_coords, disp, ii));
233  break;
234  case ROLLER:
235  rollerDataVec.push_back(
236  new RollerRigidBody(center_coords, disp, ii));
237  break;
238  case STL:
239  rollerDataVec.push_back(new STLRigidBody(center_coords, disp, ii));
240  break;
241  case NURBS:
242  // code to be executed if
243  break;
244  default:
245  // SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
246  // "this type of rigid body is not supported, check spelling "
247  // "(plane, sphere, cylinder, torus, superquadric, stl, nurbs)");
248  break;
249  }
250  CHKERR rollerDataVec.back().getBodyOptions();
251 
252  char pos_data_file[255];
253  PetscBool ctg_flag = PETSC_FALSE;
254  param = "-position_data" + to_string(ii);
255  CHKERR PetscOptionsString(param.c_str(), "", "", "", pos_data_file, 255,
256  &ctg_flag);
257  if (ctg_flag) {
258  rollerDataVec.back().positionDataParamName = string(param);
259  rollerDataVec.back().originCoords.clear();
260  rollerDataVec.back().rollerDisp.clear();
261  }
262  }
263  }
264 
265  if (with_contact && rollerDataVec.empty()) {
266  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
267  "provide at least one body data for contact (-body1 "
268  "sphere -radius r and opional -move x,y) ");
269  }
270  // }
271 
272  ierr = PetscOptionsEnd();
273  CHKERRQ(ierr);
274 
276  }
277 
280 
282  young_modulus = 1;
284  H *= young_modulus_inv;
288 
292  // cn_cont = 1; // perhaps these should be always assigned by user
293  // cn_pl = 1;
294 
296  }
297 };
298 
299 extern BlockParamData *cache;
300 extern std::map<int, BlockParamData> mat_blocks; // for use in future
301 
302 struct DataFromMove {
303  size_t comp;
306  Range ents;
307 
308  DataFromMove(size_t c, VectorDouble val, Range ents_bc)
309  : comp(c), values(val), scaledValues(val), ents(ents_bc) {}
310  DataFromMove(VectorDouble val, Range ents_bc)
311  : comp(0), values(val), scaledValues(val), ents(ents_bc) {}
312 
313  DataFromMove() = delete;
314 };
315 
317  double &scale;
318  LoadScale(double &my_scale) : scale(my_scale){};
321  nf *= scale;
323  }
324 };
325 
326 /**
327  * \brief Not used at this stage. Could be used to do some calculations,
328  * before assembly of local elements.
329  */
330 struct FePrePostProcess : public FEMethod {
331 
333  boost::ptr_vector<RigidBodyData> &rollerDataVec;
334  boost::ptr_vector<DataFromMove> &bcData;
335 
336  boost::shared_ptr<DirichletDisplacementBc> dirichletBcPtr;
337  boost::shared_ptr<MethodForForceScaling> methodsOpForMove;
338  boost::shared_ptr<MethodForForceScaling> methodsOpForRollersDisp;
339  vector<boost::shared_ptr<MethodForForceScaling>> methodsOpForRollersPosition;
340 
341  boost::shared_ptr<NestedMatrixContainer> nestedContainer;
342 
343  SmartPetscObj<Mat> SEpEp;
344  SmartPetscObj<AO> aoSEpEp;
345 
346  bool scaleMat;
348 
350  // MoFEM::Interface &mField;
352  boost::shared_ptr<DirichletDisplacementBc> dirichletBcPtr;
353  BcEntMethod(DataFromBc &data_from_dirichlet_bc,
354  boost::shared_ptr<DirichletDisplacementBc> dirichlet_bc)
355  : dataFromDirichletBc(data_from_dirichlet_bc),
356  dirichletBcPtr(dirichlet_bc) {}
357  MoFEMErrorCode preProcess() { return 0; }
358  MoFEMErrorCode postProcess() { return 0; }
361  EntityHandle v = entPtr->getEnt();
362  auto &bc_it = dataFromDirichletBc;
363  for (int i = 0; i != 3; i++) {
364  CHKERR dirichletBcPtr->calculateRotationForDof(v, bc_it);
365  if (bc_it.bc_flags[i]) {
366  if (entPtr->getEntType() == MBVERTEX) {
367  entPtr->getEntFieldData()[i] = bc_it.scaled_values[i];
368  } else {
369  entPtr->getEntFieldData()[i] = 0;
370  }
371  }
372  }
373 
375  }
376  };
377 
379  boost::ptr_vector<RigidBodyData> &roller_data,
380  boost::ptr_vector<DataFromMove> &bc_data, bool scale_mat)
381  : mField(m_field), rollerDataVec(roller_data), bcData(bc_data),
382  scaleMat(scale_mat) {}
383 
385  boost::shared_ptr<NestedMatrixContainer> nested_container) {
387  nestedContainer = nested_container;
389  }
390 
393  runOnceFlag = true;
394  CHKERR VecSetOption(ts_F, VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE);
395 
396  switch (ts_ctx) {
397  case CTX_TSSETIJACOBIAN: {
398  snes_ctx = CTX_SNESSETJACOBIAN;
399  snes_B = ts_B;
400  // CHKERR MatSetOption(ts_B, MAT_NO_OFF_PROC_ZERO_ROWS, PETSC_TRUE);
401  // CHKERR MatSetOption(ts_B, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE);
402  // FIXME: this is probably not necessary
403  // this is only for keeping constant matrices
404  if (nestedContainer) {
405  auto &nested_matrices = nestedContainer->nested_matrices;
406  if (nested_matrices(0, 0)) {
407  CHKERR MatZeroEntries(nested_matrices(0, 0));
408  CHKERR MatZeroEntries(nested_matrices(0, 1));
409  CHKERR MatZeroEntries(nested_matrices(1, 1));
410  }
411  }
412  if (SEpEp)
413  CHKERR MatZeroEntries(SEpEp);
414  break;
415  }
416  case CTX_TSSETIFUNCTION: {
417  snes_ctx = CTX_SNESSETFUNCTION;
418  snes_f = ts_F;
419  break;
420  }
421  default:
422  break;
423  }
424 
425  auto time = ts_t;
426 
427  // forceDataScaled = forceData;
428  // pressureDataScaled = pressureData;
429  if (methodsOpForMove)
430  for (auto &bdata : bcData) {
431  bdata.scaledValues = bdata.values;
433  bdata.scaledValues);
434  }
435 
436  if (time <= (*cache).rollers_stop_time) {
437 
439  for (auto &roller : (*cache).rollerDataVec) {
440  roller.BodyDispScaled = roller.rollerDisp;
442  this, methodsOpForRollersDisp, roller.BodyDispScaled);
443  }
444 
445  if (!methodsOpForRollersPosition.empty())
446 
447  for (int id = 0; id != (*cache).rollerDataVec.size(); ++id) {
448  auto &roller = (*cache).rollerDataVec[id];
449  if (!roller.positionDataParamName.empty())
451  this, methodsOpForRollersPosition[id], roller.BodyDispScaled);
452  }
453  }
454 
455  auto fix_dofs_values = [this]() {
457  auto &methodsOp = dirichletBcPtr->methodsOp;
458 
459  // for fixing dirichlet bcs
460  vector<DataFromBc> bcData;
461  CHKERR dirichletBcPtr->getRotationBcFromBlock(bcData);
462  CHKERR dirichletBcPtr->getBcDataFromSetsAndBlocks(bcData);
463 
464  for (auto &bc_it : bcData) {
465  Range all_bc_ents;
467  bc_it.scaled_values);
468  bc_it.theta = bc_it.scaled_values[0];
469  BcEntMethod method(bc_it, dirichletBcPtr);
470  for (auto &ents : bc_it.bc_ents)
471  all_bc_ents.merge(ents);
472 
473  CHKERR mField.loop_entities("U", method, &all_bc_ents);
474  }
475 
477  };
478 
479  CHKERR fix_dofs_values();
480 
481  // if (snes_ctx == CTX_SNESNONE) {
482  // if (!dofsIndices.empty()) {
483  // CHKERR VecSetValues(snes_x, dofsIndices.size(),
484  // &*dofsIndices.begin(),
485  // &*dofsValues.begin(), INSERT_VALUES);
486  // }
487  // CHKERR VecAssemblyBegin(snes_x);
488  // CHKERR VecAssemblyEnd(snes_x);
489  // }
490 
492  }
493 
496  // if (ts_ctx == CTX_TSSETIJACOBIAN)
497  // CHKERR MatView(ts_B, PETSC_VIEWER_STDOUT_SELF);
498  auto assemble = [](Mat a) {
500  if (a) {
501  CHKERR MatAssemblyBegin(a, MAT_FINAL_ASSEMBLY);
502  CHKERR MatAssemblyEnd(a, MAT_FINAL_ASSEMBLY);
503  }
505  };
506 
507  if (ts_ctx == CTX_TSSETIJACOBIAN) {
508 
509  // FIXME: this is probably not necessary
510  // this is only for keeping constant matrices
511  if (nestedContainer) {
512  auto &nested_matrices = nestedContainer->nested_matrices;
513  if (nested_matrices(0, 0)) {
514  CHKERR assemble(nested_matrices(0, 0));
515  CHKERR assemble(nested_matrices(0, 1));
516  CHKERR assemble(nested_matrices(1, 1));
517  }
518 
519  // cout << " nested_matrices(0, 0) " << endl;
520  // CHKERR MatView(nested_matrices(0, 0), PETSC_VIEWER_STDOUT_SELF);
521  // cout << " nested_matrices(0, 1) " << endl;
522  // CHKERR MatView(nested_matrices(0, 1), PETSC_VIEWER_STDOUT_SELF);
523  // cout << " nested_matrices(1, 1) " << endl;
524  // CHKERR MatView(nested_matrices(1, 1), PETSC_VIEWER_STDOUT_SELF);
525  }
526  if (SEpEp)
527  CHKERR assemble(SEpEp);
528  // if (true)
529  // CHKERR MatView(ts_B, PETSC_VIEWER_DRAW_WORLD);
530  // if (SEpEp)
531  // CHKERR MatView(SEpEp, PETSC_VIEWER_DRAW_WORLD);
532 
533  }
534 
536  }
537 };
538 
539 // }; // namespace OpContactTools
std::map< int, BlockParamData > mat_blocks
static Index< 'l', 3 > l
static Index< 'k', 3 > k
static Index< 'm', 3 > m
static Index< 'M', 3 > M
static Index< 'n', 3 > n
static Index< 'j', 3 > j
static Index< 'o', 3 > o
static Index< 'L', 3 > L
static Index< 'J', 3 > J
static Number< 2 > N2
static Index< 'p', 3 > p
BlockParamData * cache
static Index< 'N', 3 > N
static Index< 'I', 3 > I
static Number< 1 > N1
static Index< 'K', 3 > K
static Index< 'i', 3 > i
static Number< 0 > N0
@ TORUS
Definition: RigidBodies.hpp:22
@ ROLLER
Definition: RigidBodies.hpp:23
@ STL
Definition: RigidBodies.hpp:24
@ PLANE
Definition: RigidBodies.hpp:18
@ CYLINDER
Definition: RigidBodies.hpp:20
@ NURBS
Definition: RigidBodies.hpp:25
@ CONE
Definition: RigidBodies.hpp:21
@ SPHERE
Definition: RigidBodies.hpp:19
@ LASTBODY
Definition: RigidBodies.hpp:26
constexpr double a
Definition: single.cpp:4
#define MoFEMFunctionReturnHot(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:460
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:359
@ MOFEM_DATA_INCONSISTENCY
Definition: definitions.h:44
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:429
#define CHKERR
Inline error check.
Definition: definitions.h:548
#define MoFEMFunctionBeginHot
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:453
virtual MoFEMErrorCode loop_entities(const std::string field_name, EntityMethod &method, Range const *const ents=nullptr, int verb=DEFAULT_VERBOSITY)=0
Loop over field entities.
double v
phase velocity of light in medium (cm/ns)
static MoFEMErrorCodeGeneric< PetscErrorCode > ierr
Definition: Exceptions.hpp:87
PetscErrorCode MoFEMErrorCode
MoFEM/PETSc error code.
Definition: Exceptions.hpp:67
UBlasMatrix< double > MatrixDouble
Definition: Types.hpp:88
UBlasVector< double > VectorDouble
Definition: Types.hpp:79
PetscErrorCode PetscOptionsGetRealArray(PetscOptions *, const char pre[], const char name[], PetscReal dval[], PetscInt *nmax, PetscBool *set)
PetscErrorCode PetscOptionsGetEList(PetscOptions *, const char pre[], const char name[], const char *const *list, PetscInt next, PetscInt *value, PetscBool *set)
MoFEMErrorCode scaleParameters()
MoFEMErrorCode getOptionsFromCommandLine()
static Tensor4< double, 3, 3, 3, 3 > Is
virtual MoFEMErrorCode scaleParameters()=0
static vector< double > gravity
double cn_cont
double poisson
double young_modulus
static double spring_stiffness
static int closest_roller_id
static Tensor1< double, 3 > omega
double scale_constraint
double density
static Tensor2< double, 3, 3 > Omega
static RigidBodyData * closest_roller
double young_modulus_inv
int getClosestRollerID(Tensor1< T1, 3 > &t_coords, Tensor1< T2, 3 > &t_disp)
double rollers_stop_time
static boost::ptr_vector< RigidBodyData > rollerDataVec
MatrixDouble mtD
static double delta
Data from Cubit blocksets.
Definition: DirichletBC.hpp:43
DataFromMove(size_t c, VectorDouble val, Range ents_bc)
DataFromMove()=delete
VectorDouble scaledValues
DataFromMove(VectorDouble val, Range ents_bc)
VectorDouble values
BcEntMethod(DataFromBc &data_from_dirichlet_bc, boost::shared_ptr< DirichletDisplacementBc > dirichlet_bc)
MoFEMErrorCode operator()()
boost::shared_ptr< DirichletDisplacementBc > dirichletBcPtr
MoFEMErrorCode preProcess()
DataFromBc & dataFromDirichletBc
MoFEMErrorCode postProcess()
Not used at this stage. Could be used to do some calculations, before assembly of local elements.
MoFEM::Interface & mField
bool scaleMat
boost::shared_ptr< MethodForForceScaling > methodsOpForRollersDisp
SmartPetscObj< AO > aoSEpEp
boost::shared_ptr< NestedMatrixContainer > nestedContainer
MoFEMErrorCode setNestedContainer(boost::shared_ptr< NestedMatrixContainer > nested_container)
MoFEMErrorCode preProcess()
MoFEMErrorCode postProcess()
FePrePostProcess(MoFEM::Interface &m_field, boost::ptr_vector< RigidBodyData > &roller_data, boost::ptr_vector< DataFromMove > &bc_data, bool scale_mat)
boost::shared_ptr< DirichletDisplacementBc > dirichletBcPtr
bool runOnceFlag
boost::shared_ptr< MethodForForceScaling > methodsOpForMove
vector< boost::shared_ptr< MethodForForceScaling > > methodsOpForRollersPosition
boost::ptr_vector< RigidBodyData > & rollerDataVec
boost::ptr_vector< DataFromMove > & bcData
SmartPetscObj< Mat > SEpEp
LoadScale(double &my_scale)
double & scale
MoFEMErrorCode scaleNf(const FEMethod *fe, VectorDouble &nf)
Class used to scale loads, f.e. in arc-length control.
static MoFEMErrorCode applyScale(const FEMethod *fe, boost::ptr_vector< MethodForForceScaling > &methods_op, VectorDouble &nf)
Deprecated interface functions.