27 #ifndef OOMPH_DISCONTINUOUS_GALERKIN_SPACE_TIME_UNSTEADY_HEAT_MIXED_ORDER_ELEMENTS_HEADER
28 #define OOMPH_DISCONTINUOUS_GALERKIN_SPACE_TIME_UNSTEADY_HEAT_MIXED_ORDER_ELEMENTS_HEADER
32 #include <oomph-lib-config.h>
50 class SpaceTimeUnsteadyHeatEquationsBase :
public virtual FiniteElement
73 template<
unsigned SPATIAL_DIM>
134 void output(std::ostream& outfile,
const unsigned& nplot);
150 void output(FILE* file_pt,
const unsigned& nplot);
156 const unsigned& nplot,
163 std::ostream& outfile,
164 const unsigned& nplot,
202 const unsigned& nplot)
const
207 std::stringstream error_stream;
208 error_stream <<
"Space-time unsteady heat elements only store a single "
209 <<
"field so i must be 0 rather than " <<
i << std::endl;
211 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
219 for (
unsigned j = 0; j < local_loop; j++)
236 std::ofstream& file_out,
238 const unsigned& nplot,
244 std::stringstream error_stream;
245 error_stream <<
"Space-time unsteady heat elements only store a single "
246 <<
"field so i must be 0 rather than " <<
i << std::endl;
248 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
256 for (
unsigned j = 0; j < local_loop; j++)
268 for (
unsigned i = 0;
i < SPATIAL_DIM;
i++)
278 (*exact_soln_pt)(spatial_coordinates, exact_soln);
281 file_out << exact_soln[0] << std::endl;
289 std::ofstream& file_out,
291 const unsigned& nplot,
298 std::stringstream error_stream;
299 error_stream <<
"Space-time unsteady heat elements only store a single "
300 <<
"field so i must be 0 rather than " <<
i << std::endl;
302 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
310 for (
unsigned j = 0; j < local_loop; j++)
316 double interpolated_t = 0.0;
325 for (
unsigned i = 0;
i < SPATIAL_DIM;
i++)
338 (*exact_soln_pt)(interpolated_t, spatial_coordinates, exact_soln);
341 file_out << exact_soln[0] << std::endl;
359 std::stringstream error_stream;
360 error_stream <<
"These unsteady heat elements only store 1 field, \n"
361 <<
"but i is currently " <<
i << std::endl;
363 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
392 double& source)
const
404 (*Source_fct_pt)(
t, x, source);
445 unsigned n_node =
nnode();
454 DShape dpsidx(n_node, SPATIAL_DIM + 1);
458 const unsigned el_dim = this->
dim();
475 for (
unsigned j = 0; j < SPATIAL_DIM; j++)
482 for (
unsigned l = 0; l < n_node; l++)
485 for (
unsigned j = 0; j < SPATIAL_DIM; j++)
488 flux[j] +=
nodal_value(l, u_nodal_index) * dpsidx(l, j);
517 unsigned n_node =
nnode();
529 double interpolated_u = 0.0;
532 for (
unsigned l = 0; l < n_node; l++)
535 interpolated_u +=
nodal_value(l, u_nodal_index) * psi[l];
539 return interpolated_u;
565 unsigned n_node =
nnode();
577 double interpolated_u = 0.0;
580 for (
unsigned l = 0; l < n_node; l++)
583 interpolated_u +=
nodal_value(
t, l, u_nodal_index) * psi[l];
587 return interpolated_u;
611 unsigned n_node =
nnode();
620 DShape dpsidx(n_node, SPATIAL_DIM + 1);
624 const unsigned el_dim = this->
dim();
641 double interpolated_dudt = 0.0;
644 for (
unsigned l = 0; l < n_node; l++)
648 nodal_value(l, u_nodal_index) * dpsidx(l, SPATIAL_DIM);
652 return interpolated_dudt;
666 DShape& dtestdx)
const = 0;
676 DShape& dtestdx)
const = 0;
686 DShape& dpsidx)
const = 0;
692 DShape& dtestdx)
const = 0;
700 const unsigned& flag);
737 template<
unsigned SPATIAL_DIM,
unsigned NNODE_1D>
739 :
public virtual QElement<SPATIAL_DIM + 1, NNODE_1D>,
746 :
QElement<SPATIAL_DIM + 1, NNODE_1D>(),
775 void output(std::ostream& outfile,
const unsigned& n_plot)
794 void output(FILE* file_pt,
const unsigned& n_plot)
805 const unsigned& n_plot,
810 outfile, n_plot, exact_soln_pt);
818 const unsigned& n_plot,
824 outfile, n_plot, time, exact_soln_pt);
874 template<
unsigned SPATIAL_DIM,
unsigned NNODE_1D>
879 double psi_values[3][NNODE_1D];
885 OneDimLagrange::shape<NNODE_1D>(
s[0], psi_values[0]);
886 OneDimLagrange::shape<NNODE_1D>(
s[1], psi_values[1]);
889 OneDimDiscontinuousGalerkinMixedOrderBasis::shape<NNODE_1D>(
s[2],
893 for (
unsigned k = 0; k < NNODE_1D; k++)
896 for (
unsigned j = 0; j < NNODE_1D; j++)
899 for (
unsigned i = 0;
i < NNODE_1D;
i++)
902 psi[index] = psi_values[0][
i] * psi_values[1][j] * psi_values[2][k];
918 template<
unsigned SPATIAL_DIM,
unsigned NNODE_1D>
925 double psi_values[3][NNODE_1D];
926 double dpsi_values[3][NNODE_1D];
932 OneDimLagrange::shape<NNODE_1D>(
s[0], psi_values[0]);
933 OneDimLagrange::shape<NNODE_1D>(
s[1], psi_values[1]);
934 OneDimLagrange::dshape<NNODE_1D>(
s[0], dpsi_values[0]);
935 OneDimLagrange::dshape<NNODE_1D>(
s[1], dpsi_values[1]);
938 OneDimDiscontinuousGalerkinMixedOrderBasis::shape<NNODE_1D>(
s[2],
940 OneDimDiscontinuousGalerkinMixedOrderBasis::dshape<NNODE_1D>(
941 s[2], dpsi_values[2]);
944 for (
unsigned k = 0; k < NNODE_1D; k++)
947 for (
unsigned j = 0; j < NNODE_1D; j++)
950 for (
unsigned i = 0;
i < NNODE_1D;
i++)
954 dpsi_values[0][
i] * psi_values[1][j] * psi_values[2][k];
958 psi_values[0][
i] * dpsi_values[1][j] * psi_values[2][k];
962 psi_values[0][
i] * psi_values[1][j] * dpsi_values[2][k];
965 psi[index] = psi_values[0][
i] * psi_values[1][j] * psi_values[2][k];
978 template<
unsigned SPATIAL_DIM,
unsigned NNODE_1D>
985 double test_values[3][NNODE_1D];
986 double dtest_values[3][NNODE_1D];
992 OneDimLagrange::shape<NNODE_1D>(
s[0], test_values[0]);
993 OneDimLagrange::shape<NNODE_1D>(
s[1], test_values[1]);
994 OneDimLagrange::dshape<NNODE_1D>(
s[0], dtest_values[0]);
995 OneDimLagrange::dshape<NNODE_1D>(
s[1], dtest_values[1]);
998 OneDimDiscontinuousGalerkinMixedOrderTest::shape<NNODE_1D>(
s[2],
1000 OneDimDiscontinuousGalerkinMixedOrderTest::dshape<NNODE_1D>(
1001 s[2], dtest_values[2]);
1004 for (
unsigned k = 0; k < NNODE_1D; k++)
1007 for (
unsigned j = 0; j < NNODE_1D; j++)
1010 for (
unsigned i = 0;
i < NNODE_1D;
i++)
1014 dtest_values[0][
i] * test_values[1][j] * test_values[2][k];
1018 test_values[0][
i] * dtest_values[1][j] * test_values[2][k];
1022 test_values[0][
i] * test_values[1][j] * dtest_values[2][k];
1026 test_values[0][
i] * test_values[1][j] * test_values[2][k];
1042 template<
unsigned SPATIAL_DIM,
unsigned NNODE_1D>
1054 const unsigned el_dim = this->dim();
1060 std::ostringstream error_message_stream;
1063 error_message_stream <<
"Need 3D space-time elements for this to work!"
1068 OOMPH_CURRENT_FUNCTION,
1069 OOMPH_EXCEPTION_LOCATION);
1074 dshape_local_ust_heat(
s, psi, dpsidx);
1081 this->local_to_eulerian_mapping(dpsidx, inverse_jacobian);
1084 this->transform_derivatives(inverse_jacobian, dpsidx);
1091 dtest_local_ust_heat(
s, test, dtestdx);
1094 this->transform_derivatives(inverse_jacobian, dtestdx);
1107 template<
unsigned SPATIAL_DIM,
unsigned NNODE_1D>
1116 const unsigned el_dim = SPATIAL_DIM + 1;
1122 for (
unsigned i = 0;
i < el_dim;
i++)
1125 s[
i] = this->integral_pt()->knot(ipt,
i);
1129 return dshape_and_dtest_eulerian_ust_heat(
s, psi, dpsidx, test, dtestdx);
1143 template<
unsigned SPATIAL_DIM,
unsigned NNODE_1D>
1146 :
public virtual QElement<SPATIAL_DIM, NNODE_1D>
1162 template<
unsigned NNODE_1D>
1181 template<
class UNSTEADY_HEAT_ELEMENT>
1201 std::stringstream error_stream;
1205 <<
"SpaceTimeUnsteadyHeatMixedOrder elements only store a single "
1206 <<
"field so fld must be 0 rather than " << fld << std::endl;
1210 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1215 unsigned nnod = this->
nnode();
1221 for (
unsigned j = 0; j < nnod; j++)
1224 data_values[j] = std::make_pair(this->
node_pt(j), fld);
1249 std::stringstream error_stream;
1253 <<
"SpaceTimeUnsteadyHeatMixedOrder elements only store a single "
1254 <<
"field so fld must be 0 rather than " << fld << std::endl;
1258 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1287 std::stringstream error_stream;
1291 <<
"SpaceTimeUnsteadyHeatMixedOrder elements only store a single "
1292 <<
"field so fld must be 0 rather than " << fld << std::endl;
1296 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1301 unsigned n_dim = this->
dim();
1304 unsigned n_node = this->
nnode();
1310 DShape dpsidx(n_node, n_dim);
1313 DShape dtestdx(n_node, n_dim);
1319 this->dshape_and_dtest_eulerian_ust_heat(
s, psi, dpsidx, test, dtestdx);
1329 const unsigned& fld,
1337 std::stringstream error_stream;
1341 <<
"SpaceTimeUnsteadyHeatMixedOrder elements only store a single "
1342 <<
"field so fld must be 0 rather than " << fld << std::endl;
1346 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1351 unsigned u_nodal_index = this->u_index_ust_heat();
1354 unsigned n_node = this->
nnode();
1360 this->
shape(s, psi);
1363 double interpolated_u = 0.0;
1366 for (
unsigned l = 0; l < n_node; l++)
1369 interpolated_u += this->
nodal_value(t, l, u_nodal_index) * psi[l];
1373 return interpolated_u;
1385 std::stringstream error_stream;
1389 <<
"SpaceTimeUnsteadyHeatMixedOrder elements only store a single "
1390 <<
"field so fld must be 0 rather than " << fld << std::endl;
1394 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1399 return this->
nnode();
1411 std::stringstream error_stream;
1414 error_stream <<
"SpaceTimeUnsteadyHeatMixedOrder elements only store a "
1415 <<
"single field so fld must be 0 rather than " << fld
1420 error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
1425 const unsigned u_nodal_index = this->u_index_ust_heat();
1434 void output(std::ostream& outfile,
const unsigned& nplot)
1437 unsigned el_dim = this->
dim();
1449 for (
unsigned iplot = 0; iplot < num_plot_points; iplot++)
1455 for (
unsigned i = 0;
i < el_dim;
i++)
1462 outfile << this->interpolated_u_ust_heat(
s) <<
" ";
1465 outfile << this->interpolated_du_dt_ust_heat(
s) <<
" ";
1472 for (
unsigned t = 1;
t < n_prev;
t++)
1475 for (
unsigned i = 0;
i < el_dim;
i++)
1486 for (
unsigned t = 1;
t < n_prev;
t++)
1489 outfile << this->interpolated_u_ust_heat(
t,
s) <<
" ";
1493 outfile << std::endl;
1506 template<
class ELEMENT>
1519 template<
class ELEMENT>
A Class for the derivatives of shape functions The class design is essentially the same as Shape,...
TimeStepper *& time_stepper_pt()
Return the pointer to the timestepper.
unsigned ntstorage() const
Return total number of doubles stored per value to record time history of each value (one for steady ...
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional QElement.
FaceGeometry()
Constructor: Call the constructor for the appropriate lower-dimensional QElement.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...
virtual unsigned nplot_points_paraview(const unsigned &nplot) const
Return the number of actual plot points for paraview plot with parameter nplot. Broken virtual; can b...
virtual void local_coordinate_of_node(const unsigned &j, Vector< double > &s) const
Get local coordinates of node j in the element; vector sets its own size (broken virtual)
Node *& node_pt(const unsigned &n)
Return a pointer to the local node n.
double nodal_value(const unsigned &n, const unsigned &i) const
Return the i-th value stored at local node n. Produces suitably interpolated values for hanging nodes...
virtual void transform_derivatives(const DenseMatrix< double > &inverse_jacobian, DShape &dbasis) const
Convert derivative w.r.t.local coordinates to derivatives w.r.t the coordinates used to assemble the ...
virtual std::string tecplot_zone_string(const unsigned &nplot) const
Return string for tecplot zone header (when plotting nplot points in each "coordinate direction")
virtual void shape(const Vector< double > &s, Shape &psi) const =0
Calculate the geometric shape functions at local coordinate s. This function must be overloaded for e...
virtual double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s.
int nodal_local_eqn(const unsigned &n, const unsigned &i) const
Return the local equation number corresponding to the i-th value at the n-th local node.
unsigned dim() const
Return the spatial dimension of the element, i.e. the number of local coordinates required to paramet...
unsigned nnode() const
Return the number of nodes.
void(* SteadyExactSolutionFctPt)(const Vector< double > &, Vector< double > &)
Function pointer for function that computes vector-valued steady "exact solution" as .
virtual void get_s_plot(const unsigned &i, const unsigned &nplot, Vector< double > &s, const bool &shifted_to_interior=false) const
Get cector of local coordinates of plot point i (when plotting nplot points in each "coordinate direc...
virtual double local_to_eulerian_mapping(const DShape &dpsids, DenseMatrix< double > &jacobian, DenseMatrix< double > &inverse_jacobian) const
Calculate the mapping from local to Eulerian coordinates, given the derivatives of the shape function...
virtual unsigned nplot_points(const unsigned &nplot) const
Return total number of plot points (when plotting nplot points in each "coordinate direction")
virtual void write_tecplot_zone_footer(std::ostream &outfile, const unsigned &nplot) const
Add tecplot zone "footer" to output stream (when plotting nplot points in each "coordinate direction"...
void(* UnsteadyExactSolutionFctPt)(const double &, const Vector< double > &, Vector< double > &)
Function pointer for function that computes Vector-valued time-dependent function as .
static DenseMatrix< double > Dummy_matrix
Empty dense matrix used as a dummy argument to combined residual and jacobian functions in the case w...
TimeStepper *& position_time_stepper_pt()
Return a pointer to the position timestepper.
An OomphLibError object which should be thrown when an run-time error is encountered....
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
Wrapper class for projectable elements. Adds "projectability" to the underlying ELEMENT.
///////////////////////////////////////////////////////////////////// ///////////////////////////////...
int local_equation(const unsigned &fld, const unsigned &j)
Return local equation number of value j in field fld.
unsigned nfields_for_projection()
Number of fields to be projected: Just one.
double jacobian_and_shape_of_field(const unsigned &fld, const Vector< double > &s, Shape &psi)
Return Jacobian of mapping and shape functions of field fld at local coordinate s.
unsigned nhistory_values_for_coordinate_projection()
Number of positional history values (Note: count includes current value!)
ProjectableUnsteadyHeatMixedOrderSpaceTimeElement()
Constructor [this was only required explicitly from gcc 4.5.2 onwards...].
void output(std::ostream &outfile, const unsigned &nplot)
Output FE representation of soln: x,t,u or x,y,t,u at n_plot^(SPATIAL_DIM+1) plot points.
unsigned nhistory_values_for_projection(const unsigned &fld)
Number of history values to be stored for fld-th field. (Note: count includes current value!...
double get_field(const unsigned &t, const unsigned &fld, const Vector< double > &s)
Return interpolated field fld at local coordinate s, at time level t (t=0: present; t>0: history valu...
unsigned nvalue_of_field(const unsigned &fld)
Return number of values in field fld: One per node.
Vector< std::pair< Data *, unsigned > > data_values_of_field(const unsigned &fld)
Specify the values associated with field fld. The information is returned in a vector of pairs which ...
/////////////////////////////////////////////////////////////////////// /////////////////////////////...
//////////////////////////////////////////////////////////////////////// ////////////////////////////...
QUnsteadyHeatMixedOrderSpaceTimeElement()
Constructor: Call constructors for QElement and SpaceTimeUnsteadyHeatMixedOrder equations.
void dtest_local_ust_heat(const Vector< double > &s, Shape &test, DShape &dtestdx) const
Test functions & derivs. w.r.t. to local coords.
unsigned required_nvalue(const unsigned &n) const
Required number of 'values' (pinned or dofs) at node n.
void dshape_local_ust_heat(const Vector< double > &s, Shape &psi, DShape &dpsidx) const
Shape functions & derivs. w.r.t. to local coords.
double dshape_and_dtest_eulerian_at_knot_ust_heat(const unsigned &ipt, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const
Shape/test functions and derivs w.r.t. to global coords at integration point ipt; return Jacobian of ...
static const unsigned Initial_Nvalue
Static array of ints to hold number of variables at nodes: Initial_Nvalue[n].
void output(std::ostream &outfile, const unsigned &n_plot)
Output function: x,t,u or x,y,t,u at n_plot^(SPATIAL_DIM+1) plot points.
void output(FILE *file_pt)
C-style output function: x,t,u or x,y,t,u.
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function: x,t,u or x,y,t,u at n_plot^(SPATIAL_DIM+1) plot points.
void output(std::ostream &outfile)
Output function: x,t,u or x,y,t,u.
void output_fct(std::ostream &outfile, const unsigned &n_plot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output function for an exact solution: x,t,u_exact or x,y,t,u_exact at n_plot^(SPATIAL_DIM+1) plot po...
void output_fct(std::ostream &outfile, const unsigned &n_plot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Output function for a time-dependent exact solution. x,t,u_exact or x,y,t,u_exact at n_plot^(SPATIAL_...
QUnsteadyHeatMixedOrderSpaceTimeElement(const QUnsteadyHeatMixedOrderSpaceTimeElement< SPATIAL_DIM, NNODE_1D > &dummy)=delete
Broken copy constructor.
double dshape_and_dtest_eulerian_ust_heat(const Vector< double > &s, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const
Shape/test functions & derivs. w.r.t. to global coords. Return Jacobian.
void shape_ust_heat(const Vector< double > &s, Shape &psi) const
Shape functions w.r.t. to local coords.
A Class for shape functions. In simple cases, the shape functions have only one index that can be tho...
Base class so that we don't need to know the dimension just to set the source function!
void(* SpaceTimeUnsteadyHeatSourceFctPt)(const double &time, const Vector< double > &x, double &u)
Function pointer to source function fct(t,x,f(x,t)) – x is a Vector!
virtual SpaceTimeUnsteadyHeatSourceFctPt & source_fct_pt()=0
Access function: Pointer to source function.
A class for all isoparametric elements that solve the SpaceTimeUnsteadyHeatMixedOrder equations.
unsigned self_test()
Self-test: Return 0 for OK.
void disable_ALE()
Disable ALE, i.e. assert the mesh is not moving – you do this at your own risk!
void get_flux(const Vector< double > &s, Vector< double > &flux) const
Get flux: flux[i]=du/dx_i.
void fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Compute element residual Vector and element Jacobian matrix (wrapper)
SpaceTimeUnsteadyHeatSourceFctPt Source_fct_pt
Pointer to source function:
double interpolated_u_ust_heat(const unsigned &t, const Vector< double > &s) const
Return FE representation of function value u(s) at local coordinate s at previous time t (t=0: presen...
void compute_norm(double &norm)
Compute norm of FE solution.
double interpolated_u_ust_heat(const Vector< double > &s) const
Return FE representation of function value u(s) at local coordinate s.
unsigned nscalar_paraview() const
Number of scalars/fields output by this element. Reimplements broken virtual function in base class.
void scalar_value_fct_paraview(std::ofstream &file_out, const unsigned &i, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt) const
Write values of the i-th scalar field at the plot points. Needs to be implemented for each new specif...
SpaceTimeUnsteadyHeatMixedOrderEquations(const SpaceTimeUnsteadyHeatMixedOrderEquations &dummy)=delete
Broken copy constructor.
void output_element_paraview(std::ofstream &outfile, const unsigned &nplot)
C-style output FE representation of soln: x,y,u or x,y,z,u at nplot^SPATIAL_DIM plot points.
static double Default_alpha_parameter
Static default value for the Alpha parameter (thermal inertia): One for natural scaling.
double interpolated_du_dt_ust_heat(const Vector< double > &s) const
Return FE representation of function value du/dt(s) at local coordinate s.
bool ALE_is_disabled
Boolean flag to indicate if ALE formulation is disabled when time-derivatives are computed....
void scalar_value_paraview(std::ofstream &file_out, const unsigned &i, const unsigned &nplot) const
Write values of the i-th scalar field at the plot points. Needs to be implemented for each new specif...
virtual double dshape_and_dtest_eulerian_at_knot_ust_heat(const unsigned &ipt, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const =0
Shape/test functions and derivs w.r.t. to global coords at integration point ipt; return Jacobian of ...
static double Default_beta_parameter
Static default value for the Beta parameter (thermal conductivity): One for natural scaling.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Compute element residual Vector (wrapper)
double *& beta_pt()
Pointer to Beta parameter (thermal conductivity)
double * Beta_pt
Pointer to Beta parameter (thermal conductivity)
SpaceTimeUnsteadyHeatMixedOrderEquations()
Constructor: Initialises the Source_fct_pt to null and sets flag to use ALE formulation of the equati...
virtual void dshape_local_ust_heat(const Vector< double > &s, Shape &psi, DShape &dpsidx) const =0
Shape functions & derivs. w.r.t. to local coords.
const double & beta() const
Beta parameter (thermal conductivity)
void output(std::ostream &outfile)
Output with default number of plot points.
std::string scalar_name_paraview(const unsigned &i) const
Name of the i-th scalar field. Default implementation returns V1 for the first one,...
virtual void shape_ust_heat(const Vector< double > &s, Shape &psi) const =0
Shape functions w.r.t. to local coords.
void scalar_value_fct_paraview(std::ofstream &file_out, const unsigned &i, const unsigned &nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt) const
Write values of the i-th scalar field at the plot points. Needs to be implemented for each new specif...
double *& alpha_pt()
Pointer to Alpha parameter (thermal inertia)
virtual void dtest_local_ust_heat(const Vector< double > &s, Shape &test, DShape &dtestdx) const =0
Test functions & derivs. w.r.t. to local coords.
double du_dt_ust_heat(const unsigned &n) const
Calculate du/dt at the n-th local node. Uses suitably interpolated value for hanging nodes.
virtual double dshape_and_dtest_eulerian_ust_heat(const Vector< double > &s, Shape &psi, DShape &dpsidx, Shape &test, DShape &dtestdx) const =0
Shape/test functions and derivs w.r.t. to global coords at local coordinate s; return Jacobian of map...
void enable_ALE()
(Re-)enable ALE, i.e. take possible mesh motion into account when evaluating the time-derivative....
void output(FILE *file_pt)
C_style output with default number of plot points.
SpaceTimeUnsteadyHeatSourceFctPt & source_fct_pt()
Access function: Pointer to source function.
virtual void fill_in_generic_residual_contribution_ust_heat(Vector< double > &residuals, DenseMatrix< double > &jacobian, const unsigned &flag)
Compute element residual Vector only (if flag=and/or element Jacobian matrix.
virtual unsigned u_index_ust_heat() const
Return the index at which the unknown value is stored. The default value, 0, is appropriate for singl...
const double & alpha() const
Alpha parameter (thermal inertia)
virtual void get_source_ust_heat(const double &t, const unsigned &ipt, const Vector< double > &x, double &source) const
Get source term at continous time t and (Eulerian) position x. Virtual so it can be overloaded in der...
SpaceTimeUnsteadyHeatSourceFctPt source_fct_pt() const
Access function: Pointer to source function. Const version.
void compute_error(std::ostream &outfile, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt, double &error, double &norm)
Get error and norm against exact solution.
double * Alpha_pt
Pointer to Alpha parameter (thermal inertia)
void output_fct(std::ostream &outfile, const unsigned &nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output exact soln: x,y,u_exact or x,y,z,u_exact at nplot^SPATIAL_DIM plot points.
unsigned ntstorage() const
Return the number of doubles required to represent history (one for steady)
std::string string(const unsigned &i)
Return the i-th string or "" if the relevant string hasn't been defined.
//////////////////////////////////////////////////////////////////// ////////////////////////////////...