linear_elasticity_elements.cc
Go to the documentation of this file.
1// LIC// ====================================================================
2// LIC// This file forms part of oomph-lib, the object-oriented,
3// LIC// multi-physics finite-element library, available
4// LIC// at http://www.oomph-lib.org.
5// LIC//
6// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
26// Non-inline functions for elements that solve the equations of linear
27// elasticity in cartesian coordinates
28
30
31
32namespace oomph
33{
34 /// Static default value for timescale ratio (1.0 -- for natural scaling)
35 template<unsigned DIM>
37
38
39 /// ///////////////////////////////////////////////////////////////
40 /// ///////////////////////////////////////////////////////////////
41 /// ///////////////////////////////////////////////////////////////
42
43 //======================================================================
44 /// Compute the strain tensor at local coordinate s
45 //======================================================================
46 template<unsigned DIM>
48 const Vector<double>& s, DenseMatrix<double>& strain) const
49 {
50#ifdef PARANOID
51 if ((strain.ncol() != DIM) || (strain.nrow() != DIM))
52 {
53 std::ostringstream error_message;
54 error_message << "Strain matrix is " << strain.ncol() << " x "
55 << strain.nrow() << ", but dimension of the equations is "
56 << DIM << std::endl;
57 throw OomphLibError(
58 error_message.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
59 }
60
61 // Find out how many position types there are
62 unsigned n_position_type = this->nnodal_position_type();
63
64 if (n_position_type != 1)
65 {
66 throw OomphLibError("LinearElasticity is not yet implemented for more "
67 "than one position type",
68 OOMPH_CURRENT_FUNCTION,
69 OOMPH_EXCEPTION_LOCATION);
70 }
71#endif
72
73
74 // Find out how many nodes there are in the element
75 unsigned n_node = nnode();
76
77 // Find the indices at which the local velocities are stored
78 unsigned u_nodal_index[DIM];
79 for (unsigned i = 0; i < DIM; i++)
80 {
81 u_nodal_index[i] = u_index_linear_elasticity(i);
82 }
83
84 // Set up memory for the shape and derivative functions
85 Shape psi(n_node);
86 DShape dpsidx(n_node, DIM);
87
88 // Call the derivatives of the shape functions
89 (void)dshape_eulerian(s, psi, dpsidx);
90
91 // Calculate interpolated values of the derivative of global position
92 DenseMatrix<double> interpolated_dudx(DIM, DIM, 0.0);
93
94 // Loop over nodes
95 for (unsigned l = 0; l < n_node; l++)
96 {
97 // Loop over velocity components
98 for (unsigned i = 0; i < DIM; i++)
99 {
100 // Get the nodal value
101 const double u_value = this->nodal_value(l, u_nodal_index[i]);
102
103 // Loop over derivative directions
104 for (unsigned j = 0; j < DIM; j++)
105 {
106 interpolated_dudx(i, j) += u_value * dpsidx(l, j);
107 }
108 }
109 }
110
111 /// Now fill in the entries of the strain tensor
112 for (unsigned i = 0; i < DIM; i++)
113 {
114 // Do upper half of matrix
115 // Note that j must be signed here for the comparison test to work
116 // Also i must be cast to an int
117 for (int j = (DIM - 1); j >= static_cast<int>(i); j--)
118 {
119 // Off diagonal terms
120 if (static_cast<int>(i) != j)
121 {
122 strain(i, j) =
123 0.5 * (interpolated_dudx(i, j) + interpolated_dudx(j, i));
124 }
125 // Diagonal terms will including growth factor when it comes back in
126 else
127 {
128 strain(i, i) = interpolated_dudx(i, i);
129 }
130 }
131 // Matrix is symmetric so just copy lower half
132 for (int j = (i - 1); j >= 0; j--)
133 {
134 strain(i, j) = strain(j, i);
135 }
136 }
137 }
138
139
140 //======================================================================
141 /// Compute the Cauchy stress tensor at local coordinate s for
142 /// displacement formulation.
143 //======================================================================
144 template<unsigned DIM>
146 const Vector<double>& s, DenseMatrix<double>& stress) const
147 {
148#ifdef PARANOID
149 if ((stress.ncol() != DIM) || (stress.nrow() != DIM))
150 {
151 std::ostringstream error_message;
152 error_message << "Stress matrix is " << stress.ncol() << " x "
153 << stress.nrow() << ", but dimension of the equations is "
154 << DIM << std::endl;
155 throw OomphLibError(
156 error_message.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
157 }
158#endif
159
160 // Get strain
161 DenseMatrix<double> strain(DIM, DIM);
162 this->get_strain(s, strain);
163
164 // Now fill in the entries of the stress tensor without exploiting
165 // symmetry -- sorry too lazy. This fct is only used for
166 // postprocessing anyway...
167 for (unsigned i = 0; i < DIM; i++)
168 {
169 for (unsigned j = 0; j < DIM; j++)
170 {
171 stress(i, j) = 0.0;
172 for (unsigned k = 0; k < DIM; k++)
173 {
174 for (unsigned l = 0; l < DIM; l++)
175 {
176 stress(i, j) += this->E(i, j, k, l) * strain(k, l);
177 }
178 }
179 }
180 }
181 }
182
183
184 //=======================================================================
185 /// Compute the residuals for the linear elasticity equations in
186 /// cartesian coordinates. Flag indicates if we want Jacobian too.
187 //=======================================================================
188 template<unsigned DIM>
191 Vector<double>& residuals, DenseMatrix<double>& jacobian, unsigned flag)
192 {
193 // Find out how many nodes there are
194 unsigned n_node = this->nnode();
195
196#ifdef PARANOID
197 // Find out how many positional dofs there are
198 unsigned n_position_type = this->nnodal_position_type();
199
200 if (n_position_type != 1)
201 {
202 throw OomphLibError("LinearElasticity is not yet implemented for more "
203 "than one position type",
204 OOMPH_CURRENT_FUNCTION,
205 OOMPH_EXCEPTION_LOCATION);
206 }
207
208 // Throw and error if an elasticity tensor has not been set
209 if (this->Elasticity_tensor_pt == 0)
210 {
211 throw OomphLibError("No elasticity tensor set.",
212 OOMPH_CURRENT_FUNCTION,
213 OOMPH_EXCEPTION_LOCATION);
214 }
215#endif
216
217 // Find the indices at which the local displacements are stored
218 unsigned u_nodal_index[DIM];
219 for (unsigned i = 0; i < DIM; i++)
220 {
221 u_nodal_index[i] = this->u_index_linear_elasticity(i);
222 }
223
224 // Timescale ratio (non-dim density)
225 double Lambda_sq = this->lambda_sq();
226
227 // Set up memory for the shape functions
228 Shape psi(n_node);
229 DShape dpsidx(n_node, DIM);
230
231 // Set the value of Nintpt -- the number of integration points
232 unsigned n_intpt = this->integral_pt()->nweight();
233
234 // Set the vector to hold the local coordinates in the element
235 Vector<double> s(DIM);
236
237 // Integer to store the local equation number
238 int local_eqn = 0, local_unknown = 0;
239
240 // Loop over the integration points
241 for (unsigned ipt = 0; ipt < n_intpt; ipt++)
242 {
243 // Assign the values of s
244 for (unsigned i = 0; i < DIM; ++i)
245 {
246 s[i] = this->integral_pt()->knot(ipt, i);
247 }
248
249 // Get the integral weight
250 double w = this->integral_pt()->weight(ipt);
251
252 // Call the derivatives of the shape functions (and get Jacobian)
253 double J = this->dshape_eulerian_at_knot(ipt, psi, dpsidx);
254
255 // Storage for Eulerian coordinates (initialised to zero)
256 Vector<double> interpolated_x(DIM, 0.0);
257
258 // Calculate interpolated values of the derivative of global position
259 // wrt lagrangian coordinates
260 DenseMatrix<double> interpolated_dudx(DIM, DIM, 0.0);
261
262 // Setup memory for accelerations (initialised to zero)
263 Vector<double> accel(DIM, 0.0);
264
265 // Calculate displacements and derivatives and lagrangian coordinates
266 for (unsigned l = 0; l < n_node; l++)
267 {
268 // Loop over displacement components (deformed position)
269 for (unsigned i = 0; i < DIM; i++)
270 {
271 // Calculate the Lagrangian coordinates and the accelerations
272 interpolated_x[i] += this->raw_nodal_position(l, i) * psi(l);
273
274 // Only compute accelerations if inertia is switched on
275 if (this->Unsteady)
276 {
277 accel[i] += this->d2u_dt2_linear_elasticity(l, i) * psi(l);
278 }
279
280 // Get the nodal displacements
281 const double u_value = this->raw_nodal_value(l, u_nodal_index[i]);
282
283 // Loop over derivative directions
284 for (unsigned j = 0; j < DIM; j++)
285 {
286 interpolated_dudx(i, j) += u_value * dpsidx(l, j);
287 }
288 }
289 }
290
291 // Get body force at current time
292 Vector<double> b(DIM);
293 this->body_force(interpolated_x, b);
294
295 // Premultiply the weights and the Jacobian
296 double W = w * J;
297
298 //=====EQUATIONS OF LINEAR ELASTICITY ========
299
300 // Loop over the test functions, nodes of the element
301 for (unsigned l = 0; l < n_node; l++)
302 {
303 // Loop over the displacement components
304 for (unsigned a = 0; a < DIM; a++)
305 {
306 // Get the equation number
307 local_eqn = this->nodal_local_eqn(l, u_nodal_index[a]);
308 /*IF it's not a boundary condition*/
309 if (local_eqn >= 0)
310 {
311 // Acceleration and body force
312 residuals[local_eqn] += (Lambda_sq * accel[a] - b[a]) * psi(l) * W;
313
314 // Stress term
315 for (unsigned b = 0; b < DIM; b++)
316 {
317 for (unsigned c = 0; c < DIM; c++)
318 {
319 for (unsigned d = 0; d < DIM; d++)
320 {
321 // Add the stress terms to the residuals
322 residuals[local_eqn] += this->E(a, b, c, d) *
323 interpolated_dudx(c, d) *
324 dpsidx(l, b) * W;
325 }
326 }
327 }
328
329 // Jacobian entries
330 if (flag)
331 {
332 // Loop over the displacement basis functions again
333 for (unsigned l2 = 0; l2 < n_node; l2++)
334 {
335 // Loop over the displacement components again
336 for (unsigned c = 0; c < DIM; c++)
337 {
338 local_unknown = this->nodal_local_eqn(l2, u_nodal_index[c]);
339 // If it's not pinned
340 if (local_unknown >= 0)
341 {
342 // Inertial term
343 if (a == c)
344 {
345 jacobian(local_eqn, local_unknown) +=
346 Lambda_sq *
347 this->node_pt(l2)->time_stepper_pt()->weight(2, 0) *
348 psi(l) * psi(l2) * W;
349 }
350
351 for (unsigned b = 0; b < DIM; b++)
352 {
353 for (unsigned d = 0; d < DIM; d++)
354 {
355 // Add the contribution to the Jacobian matrix
356 jacobian(local_eqn, local_unknown) +=
357 this->E(a, b, c, d) * dpsidx(l2, d) * dpsidx(l, b) *
358 W;
359 }
360 }
361 } // End of if not boundary condition
362 }
363 }
364 } // End of jacobian calculation
365
366 } // End of if not boundary condition
367 } // End of loop over coordinate directions
368 } // End of loop over shape functions
369 } // End of loop over integration points
370 }
371
372
373 //=======================================================================
374 /// Output exact solution x,y,[z],u,v,[w]
375 //=======================================================================
376 template<unsigned DIM>
378 std::ostream& outfile,
379 const unsigned& nplot,
381 {
382 // Vector of local coordinates
383 Vector<double> s(DIM);
384
385 // Vector for coordintes
386 Vector<double> x(DIM);
387
388 // Tecplot header info
389 outfile << this->tecplot_zone_string(nplot);
390
391 // Exact solution Vector
392 Vector<double> exact_soln(DIM);
393
394 // Loop over plot points
395 unsigned num_plot_points = this->nplot_points(nplot);
396 for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
397 {
398 // Get local coordinates of plot point
399 this->get_s_plot(iplot, nplot, s);
400
401 // Get x position as Vector
402 this->interpolated_x(s, x);
403
404 // Get exact solution at this point
405 (*exact_soln_pt)(x, exact_soln);
406
407 // Output x,y,...,u_exact,...
408 for (unsigned i = 0; i < DIM; i++)
409 {
410 outfile << x[i] << " ";
411 }
412 for (unsigned i = 0; i < DIM; i++)
413 {
414 outfile << exact_soln[i] << " ";
415 }
416 outfile << std::endl;
417 }
418
419 // Write tecplot footer (e.g. FE connectivity lists)
420 this->write_tecplot_zone_footer(outfile, nplot);
421 }
422
423
424 //=======================================================================
425 /// Output exact solution x,y,[z],u,v,[w] (unsteady version)
426 //=======================================================================
427 template<unsigned DIM>
429 std::ostream& outfile,
430 const unsigned& nplot,
431 const double& time,
433 {
434 // Vector of local coordinates
435 Vector<double> s(DIM);
436
437 // Vector for coordintes
438 Vector<double> x(DIM);
439
440 // Tecplot header info
441 outfile << this->tecplot_zone_string(nplot);
442
443 // Exact solution Vector
444 Vector<double> exact_soln(DIM);
445
446 // Loop over plot points
447 unsigned num_plot_points = this->nplot_points(nplot);
448 for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
449 {
450 // Get local coordinates of plot point
451 this->get_s_plot(iplot, nplot, s);
452
453 // Get x position as Vector
454 this->interpolated_x(s, x);
455
456 // Get exact solution at this point
457 (*exact_soln_pt)(time, x, exact_soln);
458
459 // Output x,y,...,u_exact,...
460 for (unsigned i = 0; i < DIM; i++)
461 {
462 outfile << x[i] << " ";
463 }
464 for (unsigned i = 0; i < DIM; i++)
465 {
466 outfile << exact_soln[i] << " ";
467 }
468 outfile << std::endl;
469 }
470
471 // Write tecplot footer (e.g. FE connectivity lists)
472 this->write_tecplot_zone_footer(outfile, nplot);
473 }
474
475
476 //=======================================================================
477 /// Output: x,y,[z],u,v,[w]
478 //=======================================================================
479 template<unsigned DIM>
480 void LinearElasticityEquations<DIM>::output(std::ostream& outfile,
481 const unsigned& nplot)
482 {
483 // Set output Vector
484 Vector<double> s(DIM);
485 Vector<double> x(DIM);
486 Vector<double> u(DIM);
487
488 // Tecplot header info
489 outfile << this->tecplot_zone_string(nplot);
490
491 // Loop over plot points
492 unsigned num_plot_points = this->nplot_points(nplot);
493 for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
494 {
495 // Get local coordinates of plot point
496 this->get_s_plot(iplot, nplot, s);
497
498 // Get Eulerian coordinates and displacements
499 this->interpolated_x(s, x);
500 this->interpolated_u_linear_elasticity(s, u);
501
502 // Output the x,y,..
503 for (unsigned i = 0; i < DIM; i++)
504 {
505 outfile << x[i] << " ";
506 }
507
508 // Output u,v,..
509 for (unsigned i = 0; i < DIM; i++)
510 {
511 outfile << u[i] << " ";
512 }
513
514 outfile << std::endl;
515 }
516
517 // Write tecplot footer (e.g. FE connectivity lists)
518 this->write_tecplot_zone_footer(outfile, nplot);
519 }
520
521
522 //=======================================================================
523 /// C-style output: x,y,[z],u,v,[w]
524 //=======================================================================
525 template<unsigned DIM>
527 const unsigned& nplot)
528 {
529 // Vector of local coordinates
530 Vector<double> s(DIM);
531
532 // Tecplot header info
533 fprintf(file_pt, "%s", this->tecplot_zone_string(nplot).c_str());
534
535 // Loop over plot points
536 unsigned num_plot_points = this->nplot_points(nplot);
537 for (unsigned iplot = 0; iplot < num_plot_points; iplot++)
538 {
539 // Get local coordinates of plot point
540 this->get_s_plot(iplot, nplot, s);
541
542 // Coordinates
543 for (unsigned i = 0; i < DIM; i++)
544 {
545 fprintf(file_pt, "%g ", this->interpolated_x(s, i));
546 }
547
548 // Displacement
549 for (unsigned i = 0; i < DIM; i++)
550 {
551 fprintf(file_pt, "%g ", this->interpolated_u_linear_elasticity(s, i));
552 }
553 }
554 fprintf(file_pt, "\n");
555
556 // Write tecplot footer (e.g. FE connectivity lists)
557 this->write_tecplot_zone_footer(file_pt, nplot);
558 }
559
560
561 //======================================================================
562 /// Validate against exact velocity solution
563 /// Solution is provided via function pointer.
564 /// Plot at a given number of plot points and compute L2 error
565 /// and L2 norm of velocity solution over element.
566 //=======================================================================
567 template<unsigned DIM>
569 std::ostream& outfile,
571 double& error,
572 double& norm)
573 {
574 error = 0.0;
575 norm = 0.0;
576
577 // Vector of local coordinates
578 Vector<double> s(DIM);
579
580 // Vector for coordinates
581 Vector<double> x(DIM);
582
583 // Set the value of n_intpt
584 unsigned n_intpt = this->integral_pt()->nweight();
585
586 outfile << "ZONE" << std::endl;
587
588 // Exact solution Vector (u,v,[w])
589 Vector<double> exact_soln(DIM);
590
591 // Loop over the integration points
592 for (unsigned ipt = 0; ipt < n_intpt; ipt++)
593 {
594 // Assign values of s
595 for (unsigned i = 0; i < DIM; i++)
596 {
597 s[i] = this->integral_pt()->knot(ipt, i);
598 }
599
600 // Get the integral weight
601 double w = this->integral_pt()->weight(ipt);
602
603 // Get jacobian of mapping
604 double J = this->J_eulerian(s);
605
606 // Premultiply the weights and the Jacobian
607 double W = w * J;
608
609 // Get x position as Vector
610 this->interpolated_x(s, x);
611
612 // Get exact solution at this point
613 (*exact_soln_pt)(x, exact_soln);
614
615 // Displacement error
616 for (unsigned i = 0; i < DIM; i++)
617 {
618 norm += exact_soln[i] * exact_soln[i] * W;
619 error +=
620 (exact_soln[i] - this->interpolated_u_linear_elasticity(s, i)) *
621 (exact_soln[i] - this->interpolated_u_linear_elasticity(s, i)) * W;
622 }
623
624 // Output x,y,[z]
625 for (unsigned i = 0; i < DIM; i++)
626 {
627 outfile << x[i] << " ";
628 }
629
630 // Output u_error,v_error,[w_error]
631 for (unsigned i = 0; i < DIM; i++)
632 {
633 outfile << exact_soln[i] - this->interpolated_u_linear_elasticity(s, i)
634 << " ";
635 }
636 outfile << std::endl;
637 }
638 }
639
640 //======================================================================
641 /// Validate against exact velocity solution
642 /// Solution is provided via function pointer.
643 /// Plot at a given number of plot points and compute L2 error
644 /// and L2 norm of velocity solution over element. Unsteady version
645 //=======================================================================
646 template<unsigned DIM>
648 std::ostream& outfile,
650 const double& time,
651 double& error,
652 double& norm)
653 {
654 error = 0.0;
655 norm = 0.0;
656
657 // Vector of local coordinates
658 Vector<double> s(DIM);
659
660 // Vector for coordinates
661 Vector<double> x(DIM);
662
663 // Set the value of n_intpt
664 unsigned n_intpt = this->integral_pt()->nweight();
665
666 outfile << "ZONE" << std::endl;
667
668 // Exact solution Vector (u,v,[w])
669 Vector<double> exact_soln(DIM);
670
671 // Loop over the integration points
672 for (unsigned ipt = 0; ipt < n_intpt; ipt++)
673 {
674 // Assign values of s
675 for (unsigned i = 0; i < DIM; i++)
676 {
677 s[i] = this->integral_pt()->knot(ipt, i);
678 }
679
680 // Get the integral weight
681 double w = this->integral_pt()->weight(ipt);
682
683 // Get jacobian of mapping
684 double J = this->J_eulerian(s);
685
686 // Premultiply the weights and the Jacobian
687 double W = w * J;
688
689 // Get x position as Vector
690 this->interpolated_x(s, x);
691
692 // Get exact solution at this point
693 (*exact_soln_pt)(time, x, exact_soln);
694
695 // Displacement error
696 for (unsigned i = 0; i < DIM; i++)
697 {
698 norm += exact_soln[i] * exact_soln[i] * W;
699 error +=
700 (exact_soln[i] - this->interpolated_u_linear_elasticity(s, i)) *
701 (exact_soln[i] - this->interpolated_u_linear_elasticity(s, i)) * W;
702 }
703
704 // Output x,y,[z]
705 for (unsigned i = 0; i < DIM; i++)
706 {
707 outfile << x[i] << " ";
708 }
709
710 // Output u_error,v_error,[w_error]
711 for (unsigned i = 0; i < DIM; i++)
712 {
713 outfile << exact_soln[i] - this->interpolated_u_linear_elasticity(s, i)
714 << " ";
715 }
716 outfile << std::endl;
717 }
718 }
719
720
721 // Instantiate the required elements
723 template class LinearElasticityEquations<2>;
724
725 template class QLinearElasticityElement<3, 3>;
727 template class LinearElasticityEquations<3>;
728
729
730} // namespace oomph
static char t char * s
Definition: cfortran.h:568
cstr elem_len * i
Definition: cfortran.h:603
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
Definition: shape.h:278
unsigned long nrow() const
Return the number of rows of the matrix.
Definition: matrices.h:485
unsigned long ncol() const
Return the number of columns of the matrix.
Definition: matrices.h:491
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as .
Definition: elements.h:1759
void(* UnsteadyExactSolutionFctPt)(const double &, const Vector< double > &, Vector< double > &)
Function pointer for function that computes Vector-valued time-dependent function as .
Definition: elements.h:1765
A base class for elements that solve the equations of linear elasticity in Cartesian coordinates....
static double Default_lambda_sq_value
Static default value for timescale ratio (1.0 – for natural scaling)
void get_strain(const Vector< double > &s, DenseMatrix< double > &strain) const
Return the strain tensor.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact solution x,y,[z],u,v,[w].
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Validate against exact solution. Solution is provided via function pointer. Plot at a given number of...
void output(std::ostream &outfile)
Output: x,y,[z],u,v,[w].
void get_stress(const Vector< double > &s, DenseMatrix< double > &sigma) const
Return the Cauchy stress tensor, as calculated from the elasticity tensor at specified local coordina...
virtual void fill_in_generic_contribution_to_residuals_linear_elasticity(Vector< double > &residuals, DenseMatrix< double > &jacobian, unsigned flag)
Private helper function to compute residuals and (if requested via flag) also the Jacobian matrix.
An OomphLibError object which should be thrown when an run-time error is encountered....
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Definition: shape.h:76
//////////////////////////////////////////////////////////////////// ////////////////////////////////...