29 #ifndef OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_TRACTION_ELEMENTS_HEADER
30 #define OOMPH_FOURIER_DECOMPOSED_TIME_HARMONIC_LINEAR_ELASTICITY_TRACTION_ELEMENTS_HEADER
34 #include <oomph-lib-config.h>
48 namespace TimeHarmonicFourierDecomposedLinearElasticityTractionElementHelper
55 Vector<std::complex<double>>& load)
57 unsigned n_dim = load.size();
58 for (
unsigned i = 0;
i < n_dim;
i++)
60 load[
i] = std::complex<double>(0.0, 0.0);
75 template<
class ELEMENT>
132 ELEMENT* cast_element_pt =
dynamic_cast<ELEMENT*
>(element_pt);
135 for (
unsigned i = 0;
i < n_dim + 1;
i++)
148 ->U_index_time_harmonic_fourier_decomposed_linear_elasticity_traction
151 ->u_index_time_harmonic_fourier_decomposed_linear_elasticity(
i);
194 const unsigned&
i)
const
206 void output(std::ostream& outfile,
const unsigned& n_plot)
218 void output(FILE* file_pt,
const unsigned& n_plot)
240 template<
class ELEMENT>
241 void TimeHarmonicFourierDecomposedLinearElasticityTractionElement<
243 Vector<std::complex<double>>& traction)
245 unsigned n_dim = this->nodal_dimension();
249 interpolated_x(
s, x);
253 outer_unit_normal(
s, unit_normal);
259 get_traction(ipt, x, unit_normal, traction);
267 template<
class ELEMENT>
273 unsigned n_node = nnode();
277 unsigned n_position_type = this->nnodal_position_type();
278 if (n_position_type != 1)
280 throw OomphLibError(
"TimeHarmonicFourierDecomposedLinearElasticity is "
281 "not yet implemented for more than one position type",
282 OOMPH_CURRENT_FUNCTION,
283 OOMPH_EXCEPTION_LOCATION);
288 const unsigned n_dim = this->nodal_dimension();
291 std::vector<std::complex<unsigned>> u_nodal_index(n_dim + 1);
292 for (
unsigned i = 0;
i < n_dim + 1;
i++)
304 ->U_index_time_harmonic_fourier_decomposed_linear_elasticity_traction
315 DShape dpsids(n_node, n_dim - 1);
318 unsigned n_intpt = integral_pt()->nweight();
321 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
324 double w = integral_pt()->weight(ipt);
327 dshape_local_at_knot(ipt, psi, dpsids);
336 for (
unsigned l = 0; l < n_node; l++)
339 for (
unsigned i = 0;
i < n_dim;
i++)
342 const double x_local = nodal_position(l,
i);
343 interpolated_x[
i] += x_local * psi(l);
346 for (
unsigned j = 0; j < n_dim - 1; j++)
348 interpolated_A(j,
i) += x_local * dpsids(l, j);
355 for (
unsigned i = 0;
i < n_dim - 1;
i++)
357 for (
unsigned j = 0; j < n_dim - 1; j++)
363 for (
unsigned k = 0; k < n_dim; k++)
365 A(
i, j) += interpolated_A(
i, k) * interpolated_A(j, k);
372 outer_unit_normal(ipt, interpolated_normal);
382 Adet = A(0, 0) * A(1, 1) - A(0, 1) * A(1, 0);
386 "Wrong dimension in "
387 "TimeHarmonicFourierDecomposedLinearElasticityTractionElement",
388 "TimeHarmonicFourierDecomposedLinearElasticityTractionElement::"
389 "fill_in_contribution_to_residuals()",
390 OOMPH_EXCEPTION_LOCATION);
395 double W = w * sqrt(Adet);
399 get_traction(ipt, interpolated_x, interpolated_normal, traction);
402 for (
unsigned l = 0; l < n_node; l++)
405 for (
unsigned i = 0;
i < n_dim + 1;
i++)
408 local_eqn = this->nodal_local_eqn(l, u_nodal_index[
i].real());
413 residuals[local_eqn] -=
414 traction[
i].real() * psi(l) * interpolated_x[0] *
W;
418 local_eqn = this->nodal_local_eqn(l, u_nodal_index[
i].imag());
423 residuals[local_eqn] -=
424 traction[
i].imag() * psi(l) * interpolated_x[0] *
W;
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
FaceElements are elements that coincide with the faces of higher-dimensional "bulk" elements....
int & face_index()
Index of the face (a number that uniquely identifies the face in the element)
double zeta_nodal(const unsigned &n, const unsigned &k, const unsigned &i) const
In a FaceElement, the "global" intrinsic coordinate of the element along the boundary,...
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
A general Finite Element class.
virtual void output(std::ostream &outfile)
Output the element data — typically the values at the nodes in a format suitable for post-processing.
virtual void build_face_element(const int &face_index, FaceElement *face_element_pt)
Function for building a lower dimensional FaceElement on the specified face of the FiniteElement....
unsigned nodal_dimension() const
Return the required Eulerian dimension of the nodes in this element.
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...
A class for elements that allow the imposition of an applied traction in the equations of time-harmon...
void traction(const Vector< double > &s, Vector< std::complex< double >> &traction)
Compute traction vector at specified local coordinate Should only be used for post-processing; ignore...
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function.
TimeHarmonicFourierDecomposedLinearElasticityTractionElement(FiniteElement *const &element_pt, const int &face_index)
Constructor, which takes a "bulk" element and the value of the index and its limit.
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Fill in contribution from Jacobian.
void output(std::ostream &outfile)
Output function.
Vector< std::complex< unsigned > > U_index_time_harmonic_fourier_decomposed_linear_elasticity_traction
Index at which the i-th displacement component is stored.
void(*&)(const Vector< double > &x, const Vector< double > &n, Vector< std::complex< double >> &traction) traction_fct_pt()
Reference to the traction function pointer.
void fill_in_contribution_to_residuals_time_harmonic_fourier_decomposed_linear_elasticity_traction(Vector< double > &residuals)
Helper function that actually calculates the residuals.
double zeta_nodal(const unsigned &n, const unsigned &k, const unsigned &i) const
Specify the value of nodal zeta from the face geometry The "global" intrinsic coordinate of the eleme...
void output(FILE *file_pt)
C_style output function.
void output(std::ostream &outfile, const unsigned &n_plot)
Output function.
virtual void get_traction(const unsigned &intpt, const Vector< double > &x, const Vector< double > &n, Vector< std::complex< double >> &traction)
Get the traction vector: Pass number of integration point (dummy), Eulerian coordinate and normal vec...
void(* Traction_fct_pt)(const Vector< double > &x, const Vector< double > &n, Vector< std::complex< double >> &result)
Pointer to an imposed traction function. Arguments: Eulerian coordinate; outer unit normal; applied t...
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Return the residuals.
void Zero_traction_fct(const Vector< double > &x, const Vector< double > &N, Vector< std::complex< double >> &load)
Default load function (zero traction)
//////////////////////////////////////////////////////////////////// ////////////////////////////////...