26 #ifndef OOMPH_GEN_HH_TIME_HARMONIC_LIN_ELAST_HEADER
27 #define OOMPH_GEN_HH_TIME_HARMONIC_LIN_ELAST_HEADER
32 #include "pml_helmholtz.h"
33 #include "time_harmonic_linear_elasticity.h"
46 template<
class ELASTICITY_BULK_ELEMENT,
class HELMHOLTZ_BULK_ELEMENT>
48 :
public virtual FaceGeometry<ELASTICITY_BULK_ELEMENT>,
49 public virtual FaceElement,
50 public virtual ElementWithExternalElement
63 Vector<std::complex<unsigned>>
71 Vector<double>& residuals);
77 FiniteElement*
const& element_pt,
const int& face_index)
78 : FaceGeometry<ELASTICITY_BULK_ELEMENT>(),
84 element_pt->build_face_element(face_index,
this);
89 ELASTICITY_BULK_ELEMENT* elem_pt =
90 dynamic_cast<ELASTICITY_BULK_ELEMENT*
>(element_pt);
92 if (elem_pt->dim() == 3)
95 RefineableElement* ref_el_pt =
96 dynamic_cast<RefineableElement*
>(elem_pt);
99 if (this->has_hanging_nodes())
101 throw OomphLibError(
"This flux element will not work correctly "
102 "if nodes are hanging\n",
103 OOMPH_CURRENT_FUNCTION,
104 OOMPH_EXCEPTION_LOCATION);
114 this->set_ninteraction(1);
117 unsigned n_dim = element_pt->nodal_dimension();
120 ELASTICITY_BULK_ELEMENT* cast_element_pt =
121 dynamic_cast<ELASTICITY_BULK_ELEMENT*
>(element_pt);
124 for (
unsigned i = 0; i < n_dim; i++)
127 cast_element_pt->u_index_time_harmonic_linear_elasticity(i);
141 DenseMatrix<double>& jacobian)
147 fill_in_jacobian_from_external_interaction_by_fd(residuals, jacobian);
155 const double&
q()
const
178 void output(std::ostream& outfile,
const unsigned& n_plot)
181 unsigned n_dim = this->nodal_dimension();
184 Vector<std::complex<double>> traction(n_dim);
187 const double q_value =
q();
192 unsigned n_intpt = integral_pt()->nweight();
195 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
197 Vector<double> s_int(n_dim - 1);
198 for (
unsigned i = 0; i < n_dim - 1; i++)
200 s_int[i] = integral_pt()->knot(ipt, i);
204 Vector<double> interpolated_normal(n_dim);
205 outer_unit_normal(ipt, interpolated_normal);
208 Vector<double> zeta(1);
209 interpolated_zeta(s_int, zeta);
212 HELMHOLTZ_BULK_ELEMENT* ext_el_pt =
213 dynamic_cast<HELMHOLTZ_BULK_ELEMENT*
>(external_element_pt(0, ipt));
214 Vector<double> s_ext(external_element_local_coord(0, ipt));
215 std::complex<double> u_helmholtz =
216 ext_el_pt->interpolated_u_pml_helmholtz(s_ext);
219 ext_el_pt->interpolated_u_pml_helmholtz(s_ext);
220 traction[0] = -q_value * interpolated_normal[0] * u_helmholtz;
221 traction[1] = -q_value * interpolated_normal[1] * u_helmholtz;
223 outfile << ext_el_pt->interpolated_x(s_ext, 0) <<
" "
224 << ext_el_pt->interpolated_x(s_ext, 1) <<
" "
225 << traction[0].real() <<
" " << traction[1].real() <<
" "
226 << traction[0].imag() <<
" " << traction[1].imag() <<
" "
227 << interpolated_normal[0] <<
" " << interpolated_normal[1]
228 <<
" " << u_helmholtz.real() <<
" " << u_helmholtz.imag() <<
" "
229 << interpolated_x(s_int, 0) <<
" " << interpolated_x(s_int, 1)
231 << sqrt(pow(ext_el_pt->interpolated_x(s_ext, 0) -
232 interpolated_x(s_int, 0),
234 pow(ext_el_pt->interpolated_x(s_ext, 1) -
235 interpolated_x(s_int, 1),
237 <<
" " << zeta[0] << std::endl;
245 FaceGeometry<ELASTICITY_BULK_ELEMENT>::output(file_pt);
249 void output(FILE* file_pt,
const unsigned& n_plot)
251 FaceGeometry<ELASTICITY_BULK_ELEMENT>::output(file_pt, n_plot);
259 std::ofstream outfile;
267 std::ofstream& outfile)
270 ELASTICITY_BULK_ELEMENT* bulk_elem_pt =
271 dynamic_cast<ELASTICITY_BULK_ELEMENT*
>(this->bulk_element_pt());
274 const unsigned bulk_dim = bulk_elem_pt->dim();
275 const unsigned local_dim = this->dim();
278 Vector<double> unit_normal(bulk_dim);
281 const unsigned n_intpt = integral_pt()->nweight();
284 Vector<double> s(local_dim);
285 std::complex<double> flux(0.0, 0.0);
286 Vector<std::complex<double>> u(bulk_dim);
289 if (outfile.is_open())
296 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
299 for (
unsigned i = 0; i < local_dim; i++)
301 s[i] = integral_pt()->knot(ipt, i);
304 this->outer_unit_normal(s, unit_normal);
307 double w = integral_pt()->weight(ipt);
310 double J = J_eulerian(s);
316 Vector<double> s_bulk(local_coordinate_in_bulk(s));
319 bulk_elem_pt->interpolated_u_time_harmonic_linear_elasticity(s_bulk, u);
322 std::complex<double> u_dot_n(0.0, 0.0);
323 for (
unsigned i = 0; i < bulk_dim; i++)
325 u_dot_n += u[i] * unit_normal[i];
329 if (outfile.is_open())
331 Vector<double> x(bulk_dim);
332 interpolated_x(s, x);
333 outfile << x[0] <<
" " << x[1] <<
" " << u_dot_n.real() <<
" "
334 << u_dot_n.imag() <<
"\n";
350 std::ofstream outfile;
358 std::ofstream& outfile)
361 const unsigned local_dim = this->dim();
364 const unsigned n_intpt = integral_pt()->nweight();
367 Vector<double> s(local_dim);
368 std::complex<double> flux(0.0, 0.0);
371 if (outfile.is_open())
378 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
381 for (
unsigned i = 0; i < local_dim; i++)
383 s[i] = integral_pt()->knot(ipt, i);
387 Vector<double> s_bulk(local_coordinate_in_bulk(s));
390 HELMHOLTZ_BULK_ELEMENT* ext_el_pt =
391 dynamic_cast<HELMHOLTZ_BULK_ELEMENT*
>(external_element_pt(0, ipt));
394 unsigned nnode_ext = ext_el_pt->nnode();
397 const unsigned ext_dim = ext_el_pt->dim();
400 Shape psi_ext(nnode_ext);
401 DShape dpsi_ext_dx(nnode_ext, ext_dim);
406 (void)ext_el_pt->dshape_eulerian(s_bulk, psi_ext, dpsi_ext_dx);
409 Vector<double> unit_normal(ext_dim);
412 this->outer_unit_normal(s, unit_normal);
415 double w = integral_pt()->weight(ipt);
418 double J = J_eulerian(s);
424 Vector<std::complex<double>> interpolated_dphidx(
425 ext_dim, std::complex<double>(0.0, 0.0));
426 for (
unsigned l = 0; l < nnode_ext; l++)
429 const std::complex<double> phi_value(
430 ext_el_pt->nodal_value(l, ext_el_pt->u_index_helmholtz().real()),
431 ext_el_pt->nodal_value(l, ext_el_pt->u_index_helmholtz().imag()));
434 for (
unsigned i = 0; i < ext_dim; i++)
436 interpolated_dphidx[i] += phi_value * dpsi_ext_dx(l, i);
441 std::complex<double> dphi_dn(0.0, 0.0);
442 for (
unsigned i = 0; i < ext_dim; i++)
444 dphi_dn += interpolated_dphidx[i] * unit_normal[i];
449 double max_permitted_discrepancy = 1.0e-10;
451 pow(ext_el_pt->interpolated_x(s_bulk, 0) - interpolated_x(s, 0), 2) +
452 pow(ext_el_pt->interpolated_x(s_bulk, 1) - interpolated_x(s, 1), 2));
453 if (diff > max_permitted_discrepancy)
455 std::ostringstream error_stream;
457 <<
"Position computed by external el and face element differ by "
458 << diff <<
" which is more than the acceptable tolerance "
459 << max_permitted_discrepancy << std::endl;
460 throw OomphLibError(error_stream.str(),
461 OOMPH_CURRENT_FUNCTION,
462 OOMPH_EXCEPTION_LOCATION);
467 if (outfile.is_open())
469 Vector<double> x(ext_dim);
470 interpolated_x(s, x);
471 outfile << x[0] <<
" " << x[1] <<
" " << dphi_dn.real() <<
" "
472 << dphi_dn.imag() <<
" "
494 template<
class ELASTICITY_BULK_ELEMENT,
class HELMHOLTZ_BULK_ELEMENT>
495 double TimeHarmonicLinElastLoadedByPMLHelmholtzPressureBCElement<
496 ELASTICITY_BULK_ELEMENT,
497 HELMHOLTZ_BULK_ELEMENT>::Default_Q_Value = 1.0;
503 template<
class ELASTICITY_BULK_ELEMENT,
class HELMHOLTZ_BULK_ELEMENT>
504 void TimeHarmonicLinElastLoadedByPMLHelmholtzPressureBCElement<
505 ELASTICITY_BULK_ELEMENT,
506 HELMHOLTZ_BULK_ELEMENT>::
507 fill_in_contribution_to_residuals_helmholtz_traction(
508 Vector<double>& residuals)
511 unsigned n_node = nnode();
515 unsigned n_position_type = this->nnodal_position_type();
516 if (n_position_type != 1)
518 throw OomphLibError(
"LinearElasticity is not yet implemented for more "
519 "than one position type.",
520 OOMPH_CURRENT_FUNCTION,
521 OOMPH_EXCEPTION_LOCATION);
526 const unsigned n_dim = this->nodal_dimension();
529 std::complex<unsigned> u_nodal_index[n_dim];
530 for (
unsigned i = 0; i < n_dim; i++)
533 this->U_index_time_harmonic_linear_elasticity_helmholtz_traction[i];
541 DShape dpsids(n_node, n_dim - 1);
544 const double q_value = q();
547 Vector<std::complex<double>> traction(2);
550 unsigned n_intpt = integral_pt()->nweight();
553 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
556 double w = integral_pt()->weight(ipt);
559 dshape_local_at_knot(ipt, psi, dpsids);
562 Vector<double> interpolated_x(n_dim, 0.0);
565 DenseMatrix<double> interpolated_A(n_dim - 1, n_dim, 0.0);
568 for (
unsigned l = 0; l < n_node; l++)
571 for (
unsigned i = 0; i < n_dim; i++)
574 const double x_local = nodal_position(l, i);
575 interpolated_x[i] += x_local * psi(l);
578 for (
unsigned j = 0; j < n_dim - 1; j++)
580 interpolated_A(j, i) += x_local * dpsids(l, j);
586 DenseMatrix<double> A(n_dim - 1);
587 for (
unsigned i = 0; i < n_dim - 1; i++)
589 for (
unsigned j = 0; j < n_dim - 1; j++)
595 for (
unsigned k = 0; k < n_dim; k++)
597 A(i, j) += interpolated_A(i, k) * interpolated_A(j, k);
603 Vector<double> interpolated_normal(n_dim);
604 outer_unit_normal(ipt, interpolated_normal);
614 Adet = A(0, 0) * A(1, 1) - A(0, 1) * A(1, 0);
618 "Wrong dimension in TimeHarmonicLinElastLoadedByPressureElement",
619 OOMPH_CURRENT_FUNCTION,
620 OOMPH_EXCEPTION_LOCATION);
625 double W = w * sqrt(Adet);
628 HELMHOLTZ_BULK_ELEMENT* ext_el_pt =
629 dynamic_cast<HELMHOLTZ_BULK_ELEMENT*
>(external_element_pt(0, ipt));
630 Vector<double> s_ext(external_element_local_coord(0, ipt));
633 std::complex<double> u_helmholtz =
634 ext_el_pt->interpolated_u_pml_helmholtz(s_ext);
635 traction[0] = -q_value * interpolated_normal[0] * u_helmholtz;
636 traction[1] = -q_value * interpolated_normal[1] * u_helmholtz;
639 for (
unsigned l = 0; l < n_node; l++)
642 for (
unsigned i = 0; i < n_dim; i++)
644 local_eqn = this->nodal_local_eqn(l, u_nodal_index[i].real());
649 residuals[local_eqn] -= traction[i].real() * psi(l) * W;
652 local_eqn = this->nodal_local_eqn(l, u_nodal_index[i].imag());
657 residuals[local_eqn] -= traction[i].imag() * psi(l) * W;
680 template<
class HELMHOLTZ_BULK_ELEMENT,
class ELASTICITY_BULK_ELEMENT>
682 :
public virtual FaceGeometry<HELMHOLTZ_BULK_ELEMENT>,
683 public virtual FaceElement,
684 public virtual ElementWithExternalElement
690 FiniteElement*
const& bulk_el_pt,
const int& face_index);
712 residuals, GeneralisedElement::Dummy_matrix, 0);
719 DenseMatrix<double>& jacobian)
723 residuals, jacobian, 1);
726 fill_in_jacobian_from_external_interaction_by_fd(residuals, jacobian);
738 void output(std::ostream& outfile,
const unsigned& n_plot)
743 const unsigned n_intpt = integral_pt()->nweight();
746 Vector<double> s_int(
Dim - 1);
749 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
752 for (
unsigned i = 0; i < (
Dim - 1); i++)
754 s_int[i] = integral_pt()->knot(ipt, i);
758 Vector<double> interpolated_normal(
Dim);
759 outer_unit_normal(ipt, interpolated_normal);
761 Vector<double> zeta(1);
762 interpolated_zeta(s_int, zeta);
765 ELASTICITY_BULK_ELEMENT* ext_el_pt =
766 dynamic_cast<ELASTICITY_BULK_ELEMENT*
>(external_element_pt(0, ipt));
767 Vector<double> s_ext(external_element_local_coord(0, ipt));
768 Vector<std::complex<double>> displ(2);
769 ext_el_pt->interpolated_u_time_harmonic_linear_elasticity(s_ext, displ);
774 std::complex<double> flux = (displ[0] * interpolated_normal[0] +
775 displ[1] * interpolated_normal[1]);
778 outfile << ext_el_pt->interpolated_x(s_ext, 0) <<
" "
779 << ext_el_pt->interpolated_x(s_ext, 1) <<
" "
780 << flux.real() * interpolated_normal[0] <<
" "
781 << flux.real() * interpolated_normal[1] <<
" "
782 << flux.imag() * interpolated_normal[0] <<
" "
783 << flux.imag() * interpolated_normal[1] <<
" "
784 << interpolated_normal[0] <<
" " << interpolated_normal[1]
785 <<
" " << flux.real() <<
" " << flux.imag() <<
" "
786 << interpolated_x(s_int, 0) <<
" " << interpolated_x(s_int, 1)
788 << sqrt(pow(ext_el_pt->interpolated_x(s_ext, 0) -
789 interpolated_x(s_int, 0),
791 pow(ext_el_pt->interpolated_x(s_ext, 1) -
792 interpolated_x(s_int, 1),
794 <<
" " << zeta[0] << std::endl;
802 FaceGeometry<HELMHOLTZ_BULK_ELEMENT>::output(file_pt);
806 void output(FILE* file_pt,
const unsigned& n_plot)
808 FaceGeometry<HELMHOLTZ_BULK_ELEMENT>::output(file_pt, n_plot);
821 unsigned n_node = nnode();
827 for (
unsigned i = 0; i < n_node; i++)
833 return J_eulerian(s);
845 unsigned n_node = nnode();
848 shape_at_knot(ipt, psi);
851 for (
unsigned i = 0; i < n_node; i++)
857 return J_eulerian_at_knot(ipt);
866 Vector<double>& residuals,
867 DenseMatrix<double>& jacobian,
868 const unsigned& flag);
887 template<
class HELMHOLTZ_BULK_ELEMENT,
class ELASTICITY_BULK_ELEMENT>
889 ELASTICITY_BULK_ELEMENT>::
890 PMLHelmholtzFluxFromNormalDisplacementBCElement(
891 FiniteElement*
const& bulk_el_pt,
const int& face_index)
892 : FaceGeometry<HELMHOLTZ_BULK_ELEMENT>(), FaceElement()
897 bulk_el_pt->build_face_element(face_index,
this);
902 HELMHOLTZ_BULK_ELEMENT* elem_pt =
903 dynamic_cast<HELMHOLTZ_BULK_ELEMENT*
>(bulk_el_pt);
905 if (elem_pt->dim() == 3)
908 RefineableElement* ref_el_pt =
909 dynamic_cast<RefineableElement*
>(elem_pt);
912 if (this->has_hanging_nodes())
914 throw OomphLibError(
"This flux element will not work correctly if "
915 "nodes are hanging\n",
916 OOMPH_CURRENT_FUNCTION,
917 OOMPH_EXCEPTION_LOCATION);
927 this->set_ninteraction(1);
931 Dim = this->node_pt(0)->ndim();
948 PMLHelmholtzEquations<1>* eqn_pt =
949 dynamic_cast<PMLHelmholtzEquations<1>*
>(bulk_el_pt);
953 std::string error_string =
954 "Bulk element must inherit from PMLHelmholtzEquations.";
956 "Nodes are one dimensional, but cannot cast the bulk element to\n";
957 error_string +=
"PMLHelmholtzEquations<1>\n.";
958 error_string +=
"If you desire this functionality, you must "
959 "implement it yourself\n";
962 error_string, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
976 PMLHelmholtzEquations<2>* eqn_pt =
977 dynamic_cast<PMLHelmholtzEquations<2>*
>(bulk_el_pt);
981 std::string error_string =
982 "Bulk element must inherit from PMLHelmholtzEquations.";
984 "Nodes are two dimensional, but cannot cast the bulk element to\n";
985 error_string +=
"PMLHelmholtzEquations<2>\n.";
986 error_string +=
"If you desire this functionality, you must "
987 "implement it yourself\n";
990 error_string, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1003 PMLHelmholtzEquations<3>* eqn_pt =
1004 dynamic_cast<PMLHelmholtzEquations<3>*
>(bulk_el_pt);
1008 std::string error_string =
1009 "Bulk element must inherit from PMLHelmholtzEquations.";
1010 error_string +=
"Nodes are three dimensional, but cannot cast the "
1011 "bulk element to\n";
1012 error_string +=
"PMLHelmholtzEquations<3>\n.";
1013 error_string +=
"If you desire this functionality, you must "
1014 "implement it yourself\n";
1016 throw OomphLibError(
1017 error_string, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1029 std::ostringstream error_stream;
1030 error_stream <<
"Dimension of node is " <<
Dim
1031 <<
". It should be 1,2, or 3!" << std::endl;
1033 throw OomphLibError(
1034 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1044 template<
class HELMHOLTZ_BULK_ELEMENT,
class ELASTICITY_BULK_ELEMENT>
1046 HELMHOLTZ_BULK_ELEMENT,
1047 ELASTICITY_BULK_ELEMENT>::
1048 fill_in_generic_residual_contribution_helmholtz_flux_from_displacement(
1049 Vector<double>& residuals,
1050 DenseMatrix<double>& jacobian,
1051 const unsigned& flag)
1054 const unsigned n_node = nnode();
1057 Shape psif(n_node), testf(n_node);
1060 const unsigned n_intpt = integral_pt()->nweight();
1063 Vector<double> s(Dim - 1);
1069 const std::complex<unsigned> u_index_helmholtz =
1070 U_index_helmholtz_from_displacement;
1074 for (
unsigned ipt = 0; ipt < n_intpt; ipt++)
1077 for (
unsigned i = 0; i < (Dim - 1); i++)
1079 s[i] = integral_pt()->knot(ipt, i);
1083 double w = integral_pt()->weight(ipt);
1087 double J = shape_and_test(s, psif, testf);
1093 Vector<double> interpolated_x(Dim, 0.0);
1096 for (
unsigned l = 0; l < n_node; l++)
1099 for (
unsigned i = 0; i < Dim; i++)
1101 interpolated_x[i] += nodal_position(l, i) * psif[l];
1106 Vector<double> interpolated_normal(2);
1107 outer_unit_normal(ipt, interpolated_normal);
1110 ELASTICITY_BULK_ELEMENT* ext_el_pt =
1111 dynamic_cast<ELASTICITY_BULK_ELEMENT*
>(external_element_pt(0, ipt));
1112 Vector<double> s_ext(external_element_local_coord(0, ipt));
1113 Vector<std::complex<double>> displ(2);
1114 ext_el_pt->interpolated_u_time_harmonic_linear_elasticity(s_ext, displ);
1119 std::complex<double> flux =
1120 (displ[0] * interpolated_normal[0] + displ[1] * interpolated_normal[1]);
1125 for (
unsigned l = 0; l < n_node; l++)
1127 local_eqn = nodal_local_eqn(l, u_index_helmholtz.real());
1132 residuals[local_eqn] -= flux.real() * testf[l] * W;
1139 local_eqn = nodal_local_eqn(l, u_index_helmholtz.imag());
1144 residuals[local_eqn] -= flux.imag() * testf[l] * W;
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Add the element's contribution to its residual vector and its Jacobian matrix.
void fill_in_generic_residual_contribution_helmholtz_flux_from_displacement(Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &flag)
Add the element's contribution to its residual vector. flag=1(or 0): do (or don't) compute the contri...
PMLHelmholtzFluxFromNormalDisplacementBCElement(const PMLHelmholtzFluxFromNormalDisplacementBCElement &dummy)=delete
Broken copy constructor.
void output(FILE *file_pt)
C-style output function.
double shape_and_test(const Vector< double > &s, Shape &psi, Shape &test) const
Function to compute the shape and test functions and to return the Jacobian of mapping between local ...
unsigned Dim
The spatial dimension of the problem.
void output(std::ostream &outfile)
Output function.
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function.
std::complex< unsigned > U_index_helmholtz_from_displacement
The index at which the unknown is stored at the nodes.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Broken assignment operator.
PMLHelmholtzFluxFromNormalDisplacementBCElement(FiniteElement *const &bulk_el_pt, const int &face_index)
Constructor, takes the pointer to the "bulk" element and the face index identifying the face to which...
double shape_and_test_at_knot(const unsigned &ipt, Shape &psi, Shape &test) const
Function to compute the shape and test functions and to return the Jacobian of mapping between local ...
void output(std::ostream &outfile, const unsigned &n_plot)
Output function: flux etc at Gauss points; n_plot is ignored.
A class for elements that allow the imposition of an applied traction in the equations of time-harmon...
Vector< std::complex< unsigned > > U_index_time_harmonic_linear_elasticity_helmholtz_traction
Index at which the i-th displacement component is stored.
static double Default_Q_Value
Static default value for the ratio of stress scales used in the fluid and solid equations (default is...
void output(FILE *file_pt)
C_style output function.
double *& q_pt()
Return a pointer the ratio of stress scales used to non-dimensionalise the fluid and solid equations.
std::complex< double > global_flux_contribution_from_helmholtz()
Compute the global_flux_contribution through the helmholtz elements : we compute dphidn....
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Return the residuals.
void output(std::ostream &outfile, const unsigned &n_plot)
Output function: Plot traction etc at Gauss points nplot is ignored.
void fill_in_contribution_to_residuals_helmholtz_traction(Vector< double > &residuals)
Helper function that actually calculates the residuals.
std::complex< double > global_flux_contribution_from_solid(std::ofstream &outfile)
Compute the element's contribution to the integral of the flux over the artificial boundary....
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function.
std::complex< double > global_flux_contribution_from_helmholtz(std::ofstream &outfile)
Compute the element's contribution to the integral of the flux over the artificial boundary....
const double & q() const
Return the ratio of the stress scales used to non-dimensionalise the fluid and elasticity equations....
TimeHarmonicLinElastLoadedByPMLHelmholtzPressureBCElement(FiniteElement *const &element_pt, const int &face_index)
Constructor, which takes a "bulk" element and the value of the index and its limit.
double * Q_pt
Pointer to the ratio, , of the stress used to non-dimensionalise the fluid stresses to the stress us...
std::complex< double > global_flux_contribution_from_solid()
Compute the global_flux_contribution through the template elasticity elements : we compute u....
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Fill in contribution from Jacobian.
void output(std::ostream &outfile)
Output function.