v0.15.0
Loading...
Searching...
No Matches
SpringElement.hpp
Go to the documentation of this file.
1/** \file SpringElements.hpp
2 \brief Header file for spring element implementation
3*/
4
5
6
7#ifndef __SPRINGELEMENT_HPP__
8#define __SPRINGELEMENT_HPP__
9
12using EntData = EntitiesFieldData::EntData;
15/** \brief Set of functions declaring elements and setting operators
16 * to apply spring boundary condition
17 */
18
31
33 : public boost::enable_shared_from_this<DataAtIntegrationPtsSprings> {
34
35 boost::shared_ptr<MatrixDouble> gradDispPtr =
36 boost::make_shared<MatrixDouble>();
37 boost::shared_ptr<MatrixDouble> xAtPts = boost::make_shared<MatrixDouble>();
38 boost::shared_ptr<MatrixDouble> xInitAtPts =
39 boost::make_shared<MatrixDouble>();
40
41 boost::shared_ptr<MatrixDouble> hMat = boost::make_shared<MatrixDouble>();
42 boost::shared_ptr<MatrixDouble> FMat = boost::make_shared<MatrixDouble>();
43 boost::shared_ptr<MatrixDouble> HMat = boost::make_shared<MatrixDouble>();
44 boost::shared_ptr<MatrixDouble> invHMat =
45 boost::make_shared<MatrixDouble>();
46 boost::shared_ptr<VectorDouble> detHVec =
47 boost::make_shared<VectorDouble>();
48
49 MatrixDouble tangent1;
50 MatrixDouble tangent2;
51 MatrixDouble normalVector;
52
56
58
59 EntitiesFieldData::EntData *faceRowData;
60
61 std::map<int, BlockOptionDataSprings> mapSpring;
62 // ~DataAtIntegrationPtsSprings() {}
63 DataAtIntegrationPtsSprings(MoFEM::Interface &m_field, double scale_stiffness = 1.)
64 : mField(m_field), faceRowData(nullptr), scaleStiffness(scale_stiffness) {
65
66 ierr = setBlocks();
67 CHKERRABORT(PETSC_COMM_WORLD, ierr);
68 }
69
70 MoFEMErrorCode getParameters() {
71 MoFEMFunctionBegin; // They will be overwritten by BlockData
72 PetscOptionsBegin(PETSC_COMM_WORLD, "", "Problem", "none");
73
74 PetscOptionsEnd();
76 }
77
86
87 MoFEMErrorCode setBlocks() {
89
91 if (bit->getName().compare(0, 9, "SPRING_BC") == 0) {
92
93 const int id = bit->getMeshsetId();
94 mapSpring[id].tRis.clear();
95 CHKERR mField.get_moab().get_entities_by_dimension(
96 bit->getMeshset(), 2, mapSpring[id].tRis, true);
97
98 std::vector<double> attributes;
99 bit->getAttributes(attributes);
100 if (attributes.size() < 2) {
101 SETERRQ(PETSC_COMM_WORLD, MOFEM_ATOM_TEST_INVALID,
102 "Springs should have 2 attributes but there is %ld",
103 attributes.size());
104 }
105 mapSpring[id].iD = id;
106 mapSpring[id].springStiffnessNormal = attributes[0];
107 mapSpring[id].springStiffnessTangent = attributes[1];
108
109 // Print spring blocks after being read
110 MOFEM_LOG_C("WORLD", Sev::verbose, "\nSpring block %d\n", id);
111 MOFEM_LOG_C("WORLD", Sev::verbose, "\tNormal stiffness %3.4g\n",
112 attributes[0]);
113 MOFEM_LOG_C("WORLD", Sev::verbose, "\tTangent stiffness %3.4g\n",
114 attributes[1]);
115 }
116 }
117
119 }
120
121 private:
123 };
124
125 /**
126 * @brief RHS-operator for the spring boundary condition element
127 *
128 * Integrates virtual
129 * work of springs on displacement or spatial positions and assembles
130 * components to RHS global vector.
131 *
132 */
134
135 // vector used to store force vector for each degree of freedom
136 VectorDouble nF;
137
138 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings> commonDataPtr;
141
142 /** * @brief Integrates and assembles to global RHS vector virtual work of
143 * springs
144 *
145 * Computes virtual work of springs on spatial positions or displacements
146 * for configurational changes:
147 *
148 * \f[
149 * f_s({\mathbf x}, {\mathbf X},
150 * \delta \mathbf{x}) = \int\limits_{\partial \Omega }^{} {{\delta
151 * \mathbf{x}^T} \cdot \left[ k_{\rm n} \left( {\mathbf N} \otimes {\mathbf
152 * N} \right) + k_{\rm t} \left( {\mathbf I} - {\mathbf N} \otimes {\mathbf
153 * N} \right) \right] \left( {\mathbf x} - {\mathbf X} \right) \partial
154 * \Omega }
155 * \f]
156 *
157 * where \f$ \delta \mathbf{x} \f$ is either virtual displacement or
158 * variation of spatial positions, \f$ k_{\rm n} \f$ is the stiffness of
159 * the springs normal to the surface, \f$ k_{\rm t} \f$ is the stiffness of
160 * the springs tangential to the surface, \f$ {\mathbf N} \f$ is the normal
161 * to the surface vector based on material positions, \f$
162 * {\mathbf x} \f$ is the vector of spatial positions or displacements of
163 * the surface with springs and \f$ {\mathbf X} \f$ is the vector of
164 * material positions that is zero when displacements are considered
165 *
166 */
167
168 MoFEMErrorCode doWork(int side, EntityType type, EntData &data);
169
170 /**
171 *
172 * @param common_data_ptr Pointer to the common data for
173 * spring element
174 * @param data Variable containing data for normal and
175 * tangential stiffnesses of springs attached to the current element
176 * @param field_name String of field name for
177 * spatial positions or displacements for rows
178 */
179 OpSpringFs(boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
180 &common_data_ptr,
182 const std::string field_name)
184 field_name.c_str(), OPROW),
185 commonDataPtr(common_data_ptr), dAta(data) {
186 if (field_name.compare(0, 16, "SPATIAL_POSITION") != 0)
187 is_spatial_position = false;
188 }
189 };
190
191 /**
192 * @brief LHS-operator for the springs element
193 *
194 * Integrates Springs virtual work on spatial positions or displacements
195 * \f$ f_s \f$ derivative with respect to spatial postions or displacements
196 * on surface with springs side and assembles components of the LHS vector.
197 *
198 */
200
201 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings> commonDataPtr;
203
204 MatrixDouble locKs;
205 MatrixDouble transLocKs;
206
207 /** * @brief Integrates and assembles to global RHS vector virtual work of
208 * springs
209 *
210 * Computes virtual work of springs on spatial positions or displacements
211 * for configurational changes:
212 *
213 * \f[
214 * \textrm{D} f_s({\mathbf x}, {\mathbf X},
215 * \delta \mathbf{x})[\Delta \mathbf x] = \int\limits_{\partial \Omega }^{}
216 * {{\delta \mathbf{x}^T} \cdot \left[ k_{\rm n} \left( {\mathbf N} \otimes
217 * {\mathbf N} \right) + k_{\rm t} \left( {\mathbf I} - {\mathbf N} \otimes
218 * {\mathbf N} \right) \right] \Delta {\mathbf x} \partial \Omega } \f]
219 *
220 * where \f$ \delta \mathbf{x} \f$ is either virtual displacement or the
221 * variation of spatial positions, \f$ k_{\rm n} \f$ is the stiffness of the
222 * springs normal to the surface, \f$ k_{\rm t} \f$ is the stiffness of the
223 * springs tangential to the surface, \f$ {\mathbf N} \f$ is the normal to
224 * the surface vector based on material positions, \f$
225 * {\mathbf x} \f$ is the vector of spatial positions or displacements of
226 * the surface with springs and \f$ {\mathbf X} \f$ is the vector of
227 * material positions that is zero when displacements are considered
228 *
229 */
230 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
231 EntityType col_type, EntData &row_data,
232 EntData &col_data);
233 /**
234 *
235 * @param common_data_ptr Pointer to the common data for
236 * spring element
237 * @param data Variable containing data for normal and
238 * tangential stiffnesses of springs attached to the current element
239 * @param field_name String of field name for
240 * spatial positions or displacements for rows
241 */
242 OpSpringKs(boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
243 &common_data_ptr,
245 const std::string field_name)
247 field_name.c_str(), field_name.c_str(), OPROWCOL),
248 commonDataPtr(common_data_ptr), dAta(data) {
249 // sYmm = false;
250 }
251 };
252
255
256 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings> commonDataPtr;
258
259 MatrixDouble locKs;
260 MatrixDouble transLocKs;
261 /**
262 * @brief Compute part of the left-hand side
263 *
264 * Computes the linearisation of the material component
265 * with respect to a variation of material coordinates
266 * \f$(\Delta\mathbf{X})\f$:
267 *
268 * \f[
269 * \textrm{D} f^{\textrm{(face)}}_s(\mathbf{x}, \mathbf{X},
270 * \delta\mathbf{x})
271 * [\Delta\mathbf{X}] = -\int\limits_{\mathcal{T}_{\xi}} \,
272 * 0.5 (k_{\rm n} - k_{\rm s}) \cdot \left\{ \left[
273 * \frac{\partial\mathbf{X}}
274 * {\partial\xi} \cdot \left(\frac{\partial\Delta
275 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{x}\right)
276 * -\frac{\partial\mathbf{X}}
277 * {\partial\eta} \cdot \left(\frac{\partial\Delta
278 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{x}\right)\right] \otimes
279 * \dfrac{\mathbf{N}}{\|\mathbf{N}\|} + \dfrac{\mathbf{N}}{\|\mathbf{N}\|}
280 * \otimes \left[ \frac{\partial\mathbf{X}}
281 * {\partial\xi} \cdot \left(\frac{\partial\Delta
282 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{x}\right)
283 * -\frac{\partial\mathbf{X}}
284 * {\partial\eta} \cdot \left(\frac{\partial\Delta
285 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{x}\right)\right]
286 * - \dfrac{\mathbf{N} \otimes \mathbf{N}}{{\| \mathbf{N} \|}^3} \left[
287 * \frac{\partial\mathbf{X}}
288 * {\partial\xi} \cdot \left(\frac{\partial\Delta
289 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{x}\right)
290 * -\frac{\partial\mathbf{X}}
291 * {\partial\eta} \cdot \left(\frac{\partial\Delta
292 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{x}\right)\right] \mathbf{N}
293 * \right \}
294 * \textrm{d}\xi\textrm{d}\eta
295 * \\
296 * +\int\limits_{\mathcal{T}_{\xi}}
297 * 0.5 k_{\rm s} \left[
298 * \frac{\partial\mathbf{X}}
299 * {\partial\xi} \cdot \left(\frac{\partial\Delta
300 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{x}\right)
301 * -\frac{\partial\mathbf{X}}
302 * {\partial\eta} \cdot \left(\frac{\partial\Delta
303 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{x}\right)\right]
304 * {\mathbf N}^{\intercal} \cdot {\mathbf I} ( {\mathbf x} - {\mathbf X} )
305 * \textrm{d}\xi\textrm{d}\eta
306 * -\int\limits_{\mathcal{T}_{\xi}}
307 * {{\delta
308 * \mathbf{x}^T} \cdot \left[ k_{\rm n} \left( {\mathbf N} \otimes {\mathbf
309 * N} \right) + k_{\rm t} \left( {\mathbf I} - {\mathbf N} \otimes {\mathbf
310 * N} \right) \right] \Delta {\mathbf X} }
311 * \textrm{d}\xi\textrm{d}\eta
312 *
313 * \f]
314 *
315 */
316 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
317 EntityType col_type, EntData &row_data,
318 EntData &col_data);
319 OpSpringKs_dX(boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
320 &common_data_ptr,
322 const std::string field_name_1,
323 const std::string field_name_2)
325 field_name_1.c_str(), field_name_2.c_str(), OPROWCOL),
326 commonDataPtr(common_data_ptr), dAta(data) {
327 sYmm = false;
328 }
329 };
330
331 /**
332 * @brief Base class for LHS-operators (material) on side volumes
333 *
334 */
336
337 MatrixDouble NN;
338
339 boost::shared_ptr<DataAtIntegrationPtsSprings> dataAtSpringPts;
340
341 VectorInt rowIndices;
342 VectorInt colIndices;
343
347
350
352
353 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
354 EntityType col_type, EntData &row_data,
355 EntData &col_data);
356 virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
359 }
360 MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data);
361
362 /**
363 * @param field_name_1 String of material positions field name
364 * @param field_name_2 String of either master or material
365 * positions field name, depending of the implementation of the child class
366 * @param data_at_spring_pts Pointer to the common data for
367 * spring element
368 *
369 */
371 const string field_name_1, const string field_name_2,
372 boost::shared_ptr<DataAtIntegrationPtsSprings> data_at_spring_pts)
373 : VolOnSideUserDataOperator(field_name_1, field_name_2,
375 dataAtSpringPts(data_at_spring_pts) {
376 sYmm = false; // This will make sure to loop over all entities
377 }
378 };
379
380 /**
381 * @brief LHS-operator (material configuration) on the side volume for spring
382 * element
383 *
384 * Computes the linearisation of the material component
385 * with respect to a variation of spatial coordinates on the side volume.
386 */
389
390 /**
391 * @brief Integrates over a face contribution from a side volume
392 *
393 * Computes linearisation of the material component
394 * with respect to a variation of spatial coordinates:
395 *
396 * \f[
397 * \textrm{D} f_s(\mathbf{x}, \mathbf{X}, \delta\mathbf{X})
398 * [\Delta\mathbf{x}] = \int\limits_{\partial \Omega }^{}
399 * \mathbf{F}^{\intercal}\cdot \left( -{{\delta \mathbf{x}^T} \cdot \left[ k_{\rm n} \left( {\mathbf N} \otimes
400 * {\mathbf N} \right) + k_{\rm t} \left( {\mathbf I} - {\mathbf N} \otimes
401 * {\mathbf N} \right) \right]\right ) \Delta {\mathbf x} \partial \Omega }
402 * -\int\limits_{\mathcal{T}_{\xi}}
403 * \left\{\left[
404 * \frac{\partial\Delta\mathbf{x}}{\partial\boldsymbol{\chi}}\,\mathbf{H}^{-1}
405 * \right]^{\intercal}\cdot \left[ k_{\rm n} \left( {\mathbf N} \otimes
406 * {\mathbf N} \right) + k_{\rm t} \left( {\mathbf I} - {\mathbf N} \otimes
407 * {\mathbf N} \right) \right] \left( {\mathbf x} - {\mathbf X}
408 * \right)\right\} \cdot \delta\mathbf{X}\, \textrm{d}\xi\textrm{d}\eta \f]
409 *
410 *
411 * where \f$ \delta \mathbf{X} \f$ is variation of
412 * material positions, \f$ k_{\rm n} \f$ is the stiffness of the springs
413 * normal to the surface, \f$ k_{\rm t} \f$ is the stiffness of the springs
414 * tangential to the surface, \f$ {\mathbf N} \f$ is the normal to the
415 * surface vector based on material positions, \f$ {\mathbf x} \f$ is the
416 * vector of spatial positions or displacements of the surface with springs,
417 * \f$ {\mathbf X} \f$ is the vector of material positions,
418 * \f$\boldsymbol{\chi}\f$ are reference coordinated,
419 * \f$\mathbf{H}\f$ is the gradient of the material map.
420 *
421 *
422 */
423 MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
424
426 const string field_name_1, const string field_name_2,
427 boost::shared_ptr<DataAtIntegrationPtsSprings> data_at_spring_pts)
428 : SpringALEMaterialVolOnSideLhs(field_name_1, field_name_2,
429 data_at_spring_pts) {
430 sYmm = false; // This will make sure to loop over all entities
431 };
432 };
433
434 /**
435 * @brief Base class for LHS-operators for pressure element (material
436 * configuration)
437 *
438 * Linearisation of the material component with respect to
439 * spatial and material coordinates consists of three parts, computed
440 * by operators working on the face and on the side volume:
441 *
442 * \f[
443 * \textrm{D} \delta W^\text{material}_p(\mathbf{x}, \mathbf{X},
444 * \delta\mathbf{x})
445 * [\Delta\mathbf{x}, \Delta\mathbf{X}] = \textrm{D} \delta
446 * W^\text{(face)}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{x})
447 * [\Delta\mathbf{X}] + \textrm{D} \delta
448 * W^\text{(side volume)}_p(\mathbf{x}, \mathbf{X}, \delta\mathbf{x})
449 * [\Delta\mathbf{x}] + \textrm{D} \delta W^\text{(side volume)}_p
450 * (\mathbf{x}, \mathbf{X}, \delta\mathbf{x}) [\Delta\mathbf{X}]
451 * \f]
452 *
453 */
454
455 /**
456 * @brief LHS-operator (material configuration) on for spring element on faces
457 *
458 * This is base struct for integrating and assembling
459 * derivatives of virtual work of springs on material positions contributing
460 * to configurational changes
461 */
463
464 boost::shared_ptr<DataAtIntegrationPtsSprings> dataAtSpringPts;
465 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideFe;
466 std::string sideFeName;
467
468 MatrixDouble NN;
469 VectorInt rowIndices;
470 VectorInt colIndices;
471
475
478
480
481 virtual MoFEMErrorCode doWork(int row_side, int col_side,
482 EntityType row_type, EntityType col_type,
483 EntData &row_data, EntData &col_data) {
486 }
487
488 virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data) {
491 }
492
493 MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data);
494
496 const string field_name_1, const string field_name_2,
497 boost::shared_ptr<DataAtIntegrationPtsSprings> data_at_spring_pts,
498 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
499 std::string &side_fe_name)
500 : FaceDataOperator(field_name_1, field_name_2,
502 dataAtSpringPts(data_at_spring_pts), sideFe(side_fe),
503 sideFeName(side_fe_name) {
504 sYmm = false; // This will make sure to loop over all entities
505 }
506 };
507
508 /**
509 * @brief LHS-operator for the pressure element (material configuration)
510 *
511 * Triggers loop over operators from the side volume
512 *
513 */
515
516 /*
517 * Triggers loop over operators from the side volume
518 *
519 */
520 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
521 EntityType col_type, EntData &row_data,
522 EntData &col_data);
523
525 const string field_name_1, const string field_name_2,
526 boost::shared_ptr<DataAtIntegrationPtsSprings> data_at_pts,
527 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
528 std::string &side_fe_name)
529 : OpSpringALEMaterialLhs(field_name_1, field_name_2, data_at_pts,
530 side_fe, side_fe_name) {
531 sYmm = false; // This will make sure to loop over all entities
532 };
533 };
534
535 /**
536 * @brief LHS-operator for the pressure element (material configuration)
537 *
538 * Computes linearisation of the material component with respect to
539 * material coordinates (also triggers a loop over operators
540 * from the side volume).
541 *
542 */
544
545 /**
546 * Integrates a contribution to the left-hand side and triggers a loop
547 * over side volume operators.
548 *
549 */
550 MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type,
551 EntityType col_type, EntData &row_data,
552 EntData &col_data);
553
554 /**
555 * @brief Compute part of the left-hand side
556 *
557 * Computes the linearisation of the material component
558 * with respect to a variation of material coordinates
559 * \f$(\Delta\mathbf{X})\f$:
560 *
561 * \f[
562 * \textrm{D} f^{\textrm{(face)}}_s(\mathbf{x}, \mathbf{X},
563 * \delta\mathbf{X})
564 * [\Delta\mathbf{X}] = -\int\limits_{\mathcal{T}_{\xi}} \,
565 * 0.5 (k_{\rm n} - k_{\rm s}) \mathbf{F}^{\intercal}\cdot \left\{ \left[
566 * \frac{\partial\mathbf{X}}
567 * {\partial\xi} \cdot \left(\frac{\partial\Delta
568 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{X}\right)
569 * -\frac{\partial\mathbf{X}}
570 * {\partial\eta} \cdot \left(\frac{\partial\Delta
571 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{X}\right)\right] \otimes
572 * \dfrac{\mathbf{N}}{\|\mathbf{N}\|} + \dfrac{\mathbf{N}}{\|\mathbf{N}\|}
573 * \otimes \left[ \frac{\partial\mathbf{X}}
574 * {\partial\xi} \cdot \left(\frac{\partial\Delta
575 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{X}\right)
576 * -\frac{\partial\mathbf{X}}
577 * {\partial\eta} \cdot \left(\frac{\partial\Delta
578 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{X}\right)\right]
579 * - \dfrac{\mathbf{N} \otimes \mathbf{N}}{{\| \mathbf{N} \|}^3} \left[
580 * \frac{\partial\mathbf{X}}
581 * {\partial\xi} \cdot \left(\frac{\partial\Delta
582 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{X}\right)
583 * -\frac{\partial\mathbf{X}}
584 * {\partial\eta} \cdot \left(\frac{\partial\Delta
585 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{X}\right)\right] \mathbf{N}
586 * \right \}
587 * \textrm{d}\xi\textrm{d}\eta
588 * \\
589 * +\int\limits_{\mathcal{T}_{\xi}}
590 * 0.5 k_{\rm s} {\mathbf{F}^{\intercal}} \left[
591 * \frac{\partial\mathbf{X}}
592 * {\partial\xi} \cdot \left(\frac{\partial\Delta
593 * \mathbf{X}}{\partial\eta}\times\delta\mathbf{X}\right)
594 * -\frac{\partial\mathbf{X}}
595 * {\partial\eta} \cdot \left(\frac{\partial\Delta
596 * \mathbf{X}}{\partial\xi}\times \delta\mathbf{X}\right)\right]
597 * {\mathbf N}^{\intercal} \cdot {\mathbf I} ( {\mathbf x} - {\mathbf X} )
598 * \textrm{d}\xi\textrm{d}\eta
599 * -\int\limits_{\mathcal{T}_{\xi}}
600 * {{\delta
601 * \mathbf{X}^T} \cdot {\mathbf{F}^{\intercal}} \left[ k_{\rm n} \left(
602 * {\mathbf N} \otimes {\mathbf N} \right) + k_{\rm t} \left( {\mathbf I} -
603 * {\mathbf N} \otimes {\mathbf N} \right) \right] \Delta {\mathbf X} }
604 * \textrm{d}\xi\textrm{d}\eta
605 *
606 * \f]
607 *
608 */
609 MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
610
612 const string field_name_1, const string field_name_2,
613 boost::shared_ptr<DataAtIntegrationPtsSprings> data_at_spring_pts,
614 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
615 std::string &side_fe_name)
616 : OpSpringALEMaterialLhs(field_name_1, field_name_2, data_at_spring_pts,
617 side_fe, side_fe_name) {
618 sYmm = false; // This will make sure to loop over all entities
619 };
620 };
621
622 /**
623 * @brief LHS-operator (material configuration) on the side volume
624 *
625 * Computes the linearisation of the material component
626 * with respect to a variation of material coordinates on the side volume.
627 *
628 */
631
632 /**
633 * @brief Integrates over a face contribution from a side volume
634 *
635 * Computes linearisation of the material component
636 * with respect to a variation of material coordinates:
637 *
638 * \f[
639 * \textrm{D} f^\text{(side volume)}_s(\mathbf{x}, \mathbf{X},
640 * \delta\mathbf{X})
641 * [\Delta\mathbf{X}] = \int\limits_{\mathcal{T}_{\xi}}
642 * \left\{\left[
643 * \mathbf{h}\,\mathbf{H}^{-1}\,\frac{\partial\Delta\mathbf{X}}
644 * {\partial\boldsymbol{\chi}}\,\mathbf{H}^{-1}
645 * \right]^{\intercal}\cdot \left[ k_{\rm n} \left( {\mathbf N} \otimes
646 * {\mathbf N} \right) + k_{\rm t} \left( {\mathbf I} - {\mathbf N} \otimes
647 * {\mathbf N} \right) \right] \left( {\mathbf x} - {\mathbf X} \right)
648 * \right\} \cdot \delta\mathbf{X}\, \textrm{d}\xi\textrm{d}\eta \f]
649 *
650 * where \f$ \delta \mathbf{X} \f$ is variation of
651 * material position, \f$ k_{\rm n} \f$ is the stiffness of the springs
652 * normal to the surface, \f$ k_{\rm t} \f$ is the stiffness of the springs
653 * tangential to the surface, \f$ {\mathbf N} \f$ is the normal to the
654 * surface vector based on material positions, \f$ {\mathbf x} \f$ is the
655 * vector of spatial positions or displacements of the surface with springs,
656 * \f$ {\mathbf X} \f$ is the vector of material positions,
657 * \f$\mathbf{h}\f$ and \f$\mathbf{H}\f$ are the gradients of the
658 * spatial and material maps, respectively, and \f$\boldsymbol{\chi}\f$ are
659 * the reference coordinates.
660 *
661 *
662 */
663 MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data);
664
666 const string field_name_1, const string field_name_2,
667 boost::shared_ptr<DataAtIntegrationPtsSprings> data_at_spring_pts)
668 : SpringALEMaterialVolOnSideLhs(field_name_1, field_name_2,
669 data_at_spring_pts) {
670 sYmm = false; // This will make sure to loop over all entities
671 };
672 };
673
674 /**
675 * @brief Operator for computing deformation gradients in side volumes
676 *
677 */
679 : public VolumeElementForcesAndSourcesCoreOnContactPrismSide::
681
682 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings> commonDataPtr;
683
685
686 MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data);
687
689 const string field_name,
690 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
691 common_data_ptr,
692 bool ho_geometry = false)
693 : VolumeElementForcesAndSourcesCoreOnContactPrismSide::UserDataOperator(
695 commonDataPtr(common_data_ptr), hoGeometry(ho_geometry) {
696 doEdges = false;
697 doQuads = false;
698 doTris = false;
699 doTets = false;
700 doPrisms = false;
701 sYmm = false;
702 };
703 };
704
705 /**
706 * @brief RHS-operator for the spring boundary condition element for ALE
707 * formulation
708 *
709 * Integrates virtual
710 * work of springs on material positions involved in configurational changes
711 * and assembles components to RHS global vector.
712 *
713 */
716
717 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings> dataAtPts;
718 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> sideFe;
719 std::string sideFeName;
720 Vec F;
723
724 VectorDouble nF;
725
726 VectorInt rowIndices;
727
728 int nbRows; ///< number of dofs on rows
729 int nbIntegrationPts; ///< number of integration points
730
731 MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data);
732
733 /**
734 * @brief Integrates and assembles to global RHS vector virtual work of
735 * springs on material positions for configurational changes for ALE
736 * formulation
737 *
738 * Computes virtual work of springs on material positions for
739 * configurational changes:
740 *
741 * \f[
742 * f_s ({\mathbf x}, {\mathbf X},
743 * \delta \mathbf{X}) = -\int\limits_{\partial \Omega }^{} {{\delta
744 * \mathbf{X}^T} \cdot {\mathbf{F}^{\intercal}} \left[ k_{\rm n} \left(
745 * {\mathbf N} \otimes {\mathbf N} \right) + k_{\rm t} \left( {\mathbf I} -
746 * {\mathbf N} \otimes {\mathbf N} \right) \right] \left( {\mathbf x} -
747 * {\mathbf X} \right) \partial \Omega } \f]
748 *
749 * where \f$ \delta \mathbf{X} \f$ is
750 * material position variation, \f$ k_{\rm n} \f$ is the stiffness of the
751 * springs normal to the surface, \f$ k_{\rm t} \f$ is the stiffness of the
752 * springs tangential to the surface, \f$ {\mathbf N} \f$ is the normal to
753 * the surface vector based on material positions, \f$ {\mathbf x} \f$ is
754 * the vector of spatial positions or displacements of the surface with
755 * springs, \f$ {\mathbf X} \f$ is the vector of material positions and
756 * finally \f$\mathbf{F}\f$ is the deformation gradient
757 *
758 *
759 * \f[
760 * \mathbf{F} = \mathbf{h}(\mathbf{x})\,\mathbf{H}(\mathbf{X})^{-1} =
761 * \frac{\partial\mathbf{x}}{\partial\boldsymbol{\chi}}
762 * \frac{\partial\boldsymbol{\chi}}{\partial\mathbf{X}}
763 * \f]
764 *
765 * where \f$\mathbf{h}\f$ and \f$\mathbf{H}\f$ are the gradients of the
766 * spatial and material maps, respectively, and \f$\boldsymbol{\chi}\f$ are
767 * the reference coordinates.
768 *
769 */
770 MoFEMErrorCode iNtegrate(EntData &row_data);
771 MoFEMErrorCode aSsemble(EntData &row_data);
772
773 /**
774 *
775 * @param material_field String of field name for
776 * material positions for rows
777 * @param data_at_pts Pointer to the common data for
778 * spring element
779 * @param side_fe Pointer to the volume finite elements
780 * adjacent to the spring face element
781 * @param side_fe_name String of 3D element adjacent to
782 * the present springs elements
783 * @param data Variable containing data for normal and
784 * tangential stiffnesses of springs attached to the current element
785 *
786 */
788 const string material_field,
789 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
790 data_at_pts,
791 boost::shared_ptr<VolumeElementForcesAndSourcesCoreOnSide> side_fe,
792 std::string &side_fe_name, MetaSpringBC::BlockOptionDataSprings &data)
793 : UserDataOperator(material_field, UserDataOperator::OPROW),
794 dataAtPts(data_at_pts), sideFe(side_fe), sideFeName(side_fe_name),
795 dAta(data){};
796 };
797
798 /// \brief Computes, for material configuration, tangent vectors to face that
799 /// lie on a surface with springs
802
803 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
805 int ngp;
806 /**
807 *
808 * @param field_name String of field name for
809 * material positions for rows
810 * @param data_at_integration_pts Pointer to the common data for
811 * spring element
812 *
813 */
815 const string field_name,
816 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
820
821 /**
822 * @brief Evaluates the two tangent vector to the triangle on surface with
823 * springs based on material base coordinates
824 *
825 * Computes the two tangent vectors,\f$ {\mathbf T}^{(1)} \f$ and \f$
826 * {\mathbf T}^{(2)}\f$, based on material base coordinates based on mesh
827 * (moab vertices) coordinates:
828 *
829 * \f[
830 * {\mathbf T}^{(1)}({\mathbf X}(\xi, \eta)) =
831 * \frac{\partial\mathbf{X}(\xi,
832 * \eta)}{\partial\xi} \;\;\;\; {\mathbf T}^{(2)}({\mathbf X}(\xi, \eta)) =
833 * \frac{\partial \mathbf{X}(\xi, \eta)}
834 * {\partial\eta}
835 * \f]
836 *
837 * where \f${\mathbf X}(\xi, \eta)\f$ is the vector of material
838 * coordinates at the gauss point of surface with springs having parent
839 * coordinates \f$\xi\f$ and \f$\eta\f$ evaluated according to
840 *
841 * \f[
842 * {\mathbf X}(\xi, \eta) =
843 * \sum\limits^{3}_{i = 1}
844 * N_i(\xi, \eta){\overline{\mathbf X}}_i
845 * \f]
846 *
847 * where \f$ N_i \f$ is the shape function corresponding to the \f$
848 * i-{\rm{th}}\f$ degree of freedom in the material configuration
849 * \f${\overline{\mathbf X}}_i\f$ corresponding to the 3 nodes of the
850 * triangular face.
851 *
852 */
853
854 MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data);
855 };
856
857 /// \brief Computes, for material configuration, normal to face that lies
858 /// on a surface with springs
861
862 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
864 int ngp;
865
866 /**
867 *
868 * @param field_name String of field name for
869 * material positions for rows
870 * @param data_at_integration_pts Pointer to the common data for
871 * spring element
872 *
873 */
875 const string field_name,
876 boost::shared_ptr<MetaSpringBC::DataAtIntegrationPtsSprings>
877 data_at_integration_pts)
879 dataAtIntegrationPts(data_at_integration_pts) {}
880
881 /**
882 * @brief Evaluates normal vector to the triangle on surface with springs
883 * based on material base coordinates
884 *
885 * Computes normal vector based on material base coordinates based on mesh
886 * (moab vertices) coordinates:
887 *
888 * \f[
889 * {\mathbf N}({\mathbf X}(\xi, \eta)) =
890 * \frac{\partial\mathbf{X}(\xi,
891 * \eta)}{\partial\xi}\times\frac{\partial \mathbf{X}(\xi, \eta)}
892 * {\partial\eta}
893 * \f]
894 *
895 * where \f${\mathbf X}(\xi, \eta)\f$ is the vector of material
896 * coordinates at the gauss point of surface with springs having parent
897 * coordinates \f$\xi\f$ and \f$\eta\f$ evaluated according to
898 *
899 * \f[
900 * {\mathbf X}(\xi, \eta) =
901 * \sum\limits^{3}_{i = 1}
902 * N_i(\xi, \eta){\overline{\mathbf X}}_i
903 * \f]
904 *
905 * where \f$ N_i \f$ is the shape function corresponding to the \f$
906 * i-{\rm{th}}\f$ degree of freedom in the material configuration
907 * \f${\overline{\mathbf X}}_i\f$ corresponding to the 3 nodes of the
908 * triangular face.
909 *
910 */
911
912 MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data);
913 };
914
915 /**
916 * \brief Declare spring element
917 *
918 * Search cubit sidesets and blocksets with spring bc and declare surface
919 * element
920
921 * Blockset has to have name “SPRING_BC”. The first two attributes of the
922 * blockset are spring stiffness value.
923
924 *
925 * @param m_field Interface instance
926 * @param field_name Field name (e.g. SPATIAL_POSITION)
927 * @param mesh_nodals_positions Name of field on which ho-geometry is defined
928 * @return Error code
929 */
930 static MoFEMErrorCode addSpringElements(
931 MoFEM::Interface &m_field, const std::string field_name,
932 const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS");
933
934 /**
935 * \brief Declare spring element in ALE
936 *
937 * Search cubit sidesets and blocksets with spring bc and declare surface
938 * element
939
940 * Blockset has to have name “SPRING_BC”. The first two attributes of the
941 * blockset are spring stiffness value.
942
943 *
944 * @param m_field Interface instance
945 * @param field_name Field name (e.g. SPATIAL_POSITION)
946 * @param mesh_nodals_positions Name of field on which ho-geometry is defined
947 * @return Error code
948 */
949 static MoFEMErrorCode
950 addSpringElementsALE(MoFEM::Interface &m_field, const std::string field_name,
951 const std::string mesh_nodals_positions,
952 Range &spring_triangles);
953
954 /**
955 * \brief Implementation of spring element. Set operators to calculate LHS and
956 * RHS
957 *
958 * @param m_field Interface instance
959 * @param fe_spring_lhs_ptr Pointer to the FE instance for LHS
960 * @param fe_spring_rhs_ptr Pointer to the FE instance for RHS
961 * @param field_name Field name (e.g. SPATIAL_POSITION)
962 * @param mesh_nodals_positions Name of field on which ho-geometry is defined
963 * @return Error code
964 */
965 static MoFEMErrorCode setSpringOperators(
966 MoFEM::Interface &m_field,
967 boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr,
968 boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_rhs_ptr,
969 const std::string field_name,
970 const std::string mesh_nodals_positions = "MESH_NODE_POSITIONS",
971 double stiffness_scale = 1.);
972
973 /**
974 * \brief Implementation of spring element. Set operators to calculate LHS and
975 * RHS
976 *
977 * @param m_field Interface instance
978 * @param fe_spring_lhs_ptr Pointer to the FE instance for LHS
979 * @param fe_spring_rhs_ptr Pointer to the FE instance for RHS
980 * @param field_name Field name (e.g. SPATIAL_POSITION)
981 * @param mesh_nodals_positions Name of field on which ho-geometry is defined
982 * @return Error code
983 */
984 static MoFEMErrorCode setSpringOperatorsMaterial(
985 MoFEM::Interface &m_field,
986 boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr_dx,
987 boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_lhs_ptr_dX,
988 boost::shared_ptr<FaceElementForcesAndSourcesCore> fe_spring_rhs_ptr,
989 boost::shared_ptr<DataAtIntegrationPtsSprings> data_at_integration_pts,
990 const std::string field_name, const std::string mesh_nodals_positions,
991 std::string side_fe_name);
992};
993
994#endif //__SPRINGELEMENT_HPP__
ForcesAndSourcesCore::UserDataOperator UserDataOperator
#define MOFEM_LOG_C(channel, severity, format,...)
static PetscErrorCode ierr
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
@ BLOCKSET
@ MOFEM_ATOM_TEST_INVALID
Definition definitions.h:40
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
#define CHKERR
Inline error check.
#define _IT_CUBITMESHSETS_BY_SET_TYPE_FOR_LOOP_(MESHSET_MANAGER, CUBITBCTYPE, IT)
Iterator that loops over a specific Cubit MeshSet having a particular BC meshset in a moFEM field.
auto bit
set bit
implementation of Data Operators for Forces and Sources
Definition Common.hpp:10
constexpr auto field_name
boost::shared_ptr< MatrixDouble > xInitAtPts
boost::shared_ptr< MatrixDouble > gradDispPtr
std::map< int, BlockOptionDataSprings > mapSpring
boost::shared_ptr< MatrixDouble > xAtPts
boost::shared_ptr< MatrixDouble > invHMat
MoFEMErrorCode getBlockData(BlockOptionDataSprings &data)
boost::shared_ptr< MatrixDouble > hMat
boost::shared_ptr< MatrixDouble > FMat
DataAtIntegrationPtsSprings(MoFEM::Interface &m_field, double scale_stiffness=1.)
boost::shared_ptr< MatrixDouble > HMat
boost::shared_ptr< VectorDouble > detHVec
EntitiesFieldData::EntData * faceRowData
Operator for computing deformation gradients in side volumes.
MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data)
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > commonDataPtr
OpCalculateDeformation(const string field_name, boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > common_data_ptr, bool ho_geometry=false)
Computes, for material configuration, normal to face that lies on a surface with springs.
OpGetNormalSpEle(const string field_name, boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > data_at_integration_pts)
MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data)
Evaluates normal vector to the triangle on surface with springs based on material base coordinates.
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > dataAtIntegrationPts
Computes, for material configuration, tangent vectors to face that lie on a surface with springs.
OpGetTangentSpEle(const string field_name, boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > dataAtIntegrationPts)
MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data)
Evaluates the two tangent vector to the triangle on surface with springs based on material base coord...
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > dataAtIntegrationPts
LHS-operator for the pressure element (material configuration)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
OpSpringALEMaterialLhs_dX_dX(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_spring_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Compute part of the left-hand side.
LHS-operator for the pressure element (material configuration)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
OpSpringALEMaterialLhs_dX_dx(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name)
Base class for LHS-operators for pressure element (material configuration)
boost::shared_ptr< DataAtIntegrationPtsSprings > dataAtSpringPts
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideFe
OpSpringALEMaterialLhs(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_spring_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name)
virtual MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data)
virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
RHS-operator for the spring boundary condition element for ALE formulation.
int nbIntegrationPts
number of integration points
MoFEMErrorCode aSsemble(EntData &row_data)
MoFEMErrorCode doWork(int side, EntityType type, EntData &row_data)
Operator for linear form, usually to calculate values on right hand side.
MetaSpringBC::BlockOptionDataSprings & dAta
boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > sideFe
int nbRows
number of dofs on rows
MoFEMErrorCode iNtegrate(EntData &row_data)
Integrates and assembles to global RHS vector virtual work of springs on material positions for confi...
OpSpringFsMaterial(const string material_field, boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > data_at_pts, boost::shared_ptr< VolumeElementForcesAndSourcesCoreOnSide > side_fe, std::string &side_fe_name, MetaSpringBC::BlockOptionDataSprings &data)
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > dataAtPts
RHS-operator for the spring boundary condition element.
MetaSpringBC::BlockOptionDataSprings & dAta
OpSpringFs(boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > &common_data_ptr, MetaSpringBC::BlockOptionDataSprings &data, const std::string field_name)
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > commonDataPtr
MoFEMErrorCode doWork(int side, EntityType type, EntData &data)
Integrates and assembles to global RHS vector virtual work of springs.
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Compute part of the left-hand side.
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > commonDataPtr
MetaSpringBC::BlockOptionDataSprings & dAta
OpSpringKs_dX(boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > &common_data_ptr, MetaSpringBC::BlockOptionDataSprings &data, const std::string field_name_1, const std::string field_name_2)
LHS-operator for the springs element.
MetaSpringBC::BlockOptionDataSprings & dAta
OpSpringKs(boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > &common_data_ptr, MetaSpringBC::BlockOptionDataSprings &data, const std::string field_name)
boost::shared_ptr< MetaSpringBC::DataAtIntegrationPtsSprings > commonDataPtr
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Integrates and assembles to global RHS vector virtual work of springs.
LHS-operator (material configuration) on the side volume.
SpringALEMaterialVolOnSideLhs_dX_dX(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_spring_pts)
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Integrates over a face contribution from a side volume.
LHS-operator (material configuration) on the side volume for spring element.
MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
Integrates over a face contribution from a side volume.
SpringALEMaterialVolOnSideLhs_dX_dx(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_spring_pts)
Base class for LHS-operators (material) on side volumes.
SpringALEMaterialVolOnSideLhs(const string field_name_1, const string field_name_2, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_spring_pts)
virtual MoFEMErrorCode iNtegrate(EntData &row_data, EntData &col_data)
boost::shared_ptr< DataAtIntegrationPtsSprings > dataAtSpringPts
MoFEMErrorCode aSsemble(EntData &row_data, EntData &col_data)
MoFEMErrorCode doWork(int row_side, int col_side, EntityType row_type, EntityType col_type, EntData &row_data, EntData &col_data)
Operator for bi-linear form, usually to calculate values on left hand side.
Set of functions declaring elements and setting operators to apply spring boundary condition.
static MoFEMErrorCode setSpringOperators(MoFEM::Interface &m_field, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_rhs_ptr, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS", double stiffness_scale=1.)
Implementation of spring element. Set operators to calculate LHS and RHS.
static MoFEMErrorCode addSpringElementsALE(MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions, Range &spring_triangles)
Declare spring element in ALE.
static MoFEMErrorCode addSpringElements(MoFEM::Interface &m_field, const std::string field_name, const std::string mesh_nodals_positions="MESH_NODE_POSITIONS")
Declare spring element.
static MoFEMErrorCode setSpringOperatorsMaterial(MoFEM::Interface &m_field, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr_dx, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_lhs_ptr_dX, boost::shared_ptr< FaceElementForcesAndSourcesCore > fe_spring_rhs_ptr, boost::shared_ptr< DataAtIntegrationPtsSprings > data_at_integration_pts, const std::string field_name, const std::string mesh_nodals_positions, std::string side_fe_name)
Implementation of spring element. Set operators to calculate LHS and RHS.
virtual moab::Interface & get_moab()=0
bool sYmm
If true assume that matrix is symmetric structure.
Deprecated interface functions.
Data on single entity (This is passed as argument to DataOperator::doWork)
@ OPCOL
operator doWork function is executed on FE columns
@ OPROW
operator doWork function is executed on FE rows
@ OPROWCOL
operator doWork is executed on FE rows &columns