26 #ifndef AXISYMMETRIC_ADVECTION_NAVIER_STOKES_ELEMENTS_HEADER
27 #define AXISYMMETRIC_ADVECTION_NAVIER_STOKES_ELEMENTS_HEADER
32 #include "axisym_advection_diffusion.h"
33 #include "axisym_navier_stokes.h"
58 public virtual QAxisymAdvectionDiffusionElement<3>,
59 public virtual AxisymmetricQCrouzeixRaviartElement
68 QAxisymAdvectionDiffusionElement<3>(),
69 AxisymmetricQCrouzeixRaviartElement() {}
76 {
return (QAxisymAdvectionDiffusionElement<3>::required_nvalue(n) +
77 AxisymmetricQCrouzeixRaviartElement::required_nvalue(n));}
83 AxisymmetricNavierStokesEquations::disable_ALE();
84 AxisymAdvectionDiffusionEquations::enable_ALE();
91 AxisymmetricNavierStokesEquations::enable_ALE();
92 AxisymAdvectionDiffusionEquations::enable_ALE();
100 return AxisymmetricNavierStokesEquations::nscalar_paraview()
101 + AxisymAdvectionDiffusionEquations::nscalar_paraview();
108 const unsigned& nplot)
const
110 AxisymmetricNavierStokesEquations::scalar_value_paraview(file_out,i,nplot);
111 AxisymAdvectionDiffusionEquations::scalar_value_paraview(file_out,i,nplot);
117 return AxisymmetricNavierStokesEquations::scalar_name_paraview(i);
124 void output(std::ostream &outfile) {FiniteElement::output(outfile);}
129 void output(std::ostream &outfile,
const unsigned &nplot)
135 outfile << this->tecplot_zone_string(nplot);
138 unsigned num_plot_points=this->nplot_points(nplot);
139 for (
unsigned iplot=0;iplot<num_plot_points;iplot++)
142 this->get_s_plot(iplot,nplot,s);
145 for(
unsigned i=0;i<2;i++)
146 {outfile << this->interpolated_x(s,i) <<
" ";}
149 for(
unsigned i=0;i<3;i++)
150 {outfile << this->interpolated_u_axi_nst(s,i) <<
" ";}
153 outfile << this->interpolated_p_axi_nst(s) <<
" ";
156 outfile << this->interpolated_u_axi_adv_diff(s) <<
" ";
160 outfile << std::endl;
163 this->write_tecplot_zone_footer(outfile,nplot);
169 {FiniteElement::output(file_pt);}
172 void output(FILE* file_pt,
const unsigned &n_plot)
173 {FiniteElement::output(file_pt,n_plot);}
176 void output_fct(std::ostream &outfile,
const unsigned &Nplot,
177 FiniteElement::SteadyExactSolutionFctPt
179 {FiniteElement::output_fct(outfile,Nplot,exact_soln_pt);}
184 void output_fct(std::ostream &outfile,
const unsigned &Nplot,
186 FiniteElement::UnsteadyExactSolutionFctPt
190 output_fct(outfile,Nplot,time,exact_soln_pt);
201 {AxisymmetricQCrouzeixRaviartElement::compute_norm(norm);}
210 FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt,
212 double& error,
double& norm)
213 {FiniteElement::compute_error(outfile,exact_soln_pt,
222 FiniteElement::SteadyExactSolutionFctPt exact_soln_pt,
223 double& error,
double& norm)
224 {FiniteElement::compute_error(outfile,exact_soln_pt,error,norm);}
230 const Vector<double> &s,
const Vector<double>& x,
231 Vector<double>& wind)
const
234 this->interpolated_u_axi_nst(s,wind);
245 AxisymmetricNavierStokesEquations::fill_in_contribution_to_residuals(residuals);
248 AxisymAdvectionDiffusionEquations::fill_in_contribution_to_residuals(residuals);
252 #ifdef USE_FD_JACOBIAN_FOR_ADVECTION_DIFFUSION_NAVIER_STOKES_ELEMENT
258 DenseMatrix<double> &jacobian)
261 FiniteElement::fill_in_contribution_to_jacobian(residuals,jacobian);
269 DenseMatrix<double> &jacobian)
273 unsigned u_nodal_nst[3];
274 for(
unsigned i=0;i<3;i++)
275 {u_nodal_nst[i] = this->u_index_nst(i);}
283 unsigned n_dof = this->ndof();
286 Vector<double> newres(n_dof);
289 int local_unknown =0, local_eqn = 0;
292 double fd_step = FiniteElement::Default_fd_jacobian_step;
295 unsigned n_node = this->nnode();
301 for(
unsigned n=0;n<n_node;n++)
304 for(
unsigned i=0;i<3;i++)
307 local_unknown = this->nodal_local_eqn(n,u_nodal_nst[i]);
310 if(local_unknown >= 0)
313 double *value_pt = this->node_pt(n)->value_pt(u_nodal_nst[i]);
316 double old_var = *value_pt;
319 *value_pt += fd_step;
325 for(
unsigned m=0;m<n_dof;m++) {newres[m] = 0.0;}
326 AxisymAdvectionDiffusionEquations::
327 fill_in_contribution_to_residuals(newres);
331 for(
unsigned m=0;m<n_node;m++)
334 local_eqn = this->nodal_local_eqn(m,C_nodal_adv_diff);
339 double sum = (newres[local_eqn] - residuals[local_eqn])/fd_step;
340 jacobian(local_eqn,local_unknown) = sum;
356 DenseMatrix<double> &jacobian)
360 AxisymmetricNavierStokesEquations::
361 fill_in_contribution_to_jacobian(residuals,jacobian);
365 AxisymAdvectionDiffusionEquations::
366 fill_in_contribution_to_jacobian(residuals,jacobian);
378 Vector<double> &residuals, DenseMatrix<double> &jacobian,
379 DenseMatrix<double> &mass_matrix)
381 AxisymmetricNavierStokesEquations::
382 fill_in_contribution_to_jacobian_and_mass_matrix(
383 residuals,jacobian,mass_matrix);
385 AxisymAdvectionDiffusionEquations::
386 fill_in_contribution_to_jacobian_and_mass_matrix(
387 residuals,jacobian,mass_matrix);
401 public virtual QElement<1,3>
412 public virtual PointElement
A class that solves the Boussinesq approximation of the Navier–Stokes and energy equations by couplin...
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...
void output(std::ostream &outfile)
Overload the standard output function with the broken default.
void output_fct(std::ostream &outfile, const unsigned &Nplot, FiniteElement::SteadyExactSolutionFctPt exact_soln_pt)
Output function for an exact solution: Broken default.
void output_fct(std::ostream &outfile, const unsigned &Nplot, const double &time, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt)
Output function for a time-dependent exact solution: Broken default.
void compute_error(std::ostream &outfile, FiniteElement::UnsteadyExactSolutionFctPt exact_soln_pt, const double &time, double &error, double &norm)
Validate against exact solution at given time Solution is provided via function pointer....
unsigned required_nvalue(const unsigned &n) const
The required number of values stored at the nodes is the sum of the required values of the two single...
void output(std::ostream &outfile, const unsigned &nplot)
Output function: Output x, y, u, v, p, theta at Nplot^DIM plot points.
void fill_in_contribution_to_jacobian_and_mass_matrix(Vector< double > &residuals, DenseMatrix< double > &jacobian, DenseMatrix< double > &mass_matrix)
Add the element's contribution to its residuals vector, jacobian matrix and mass matrix.
void fill_in_contribution_to_residuals(Vector< double > &residuals)
Calculate the element's contribution to the residual vector. Recall that fill_in_* functions MUST NOT...
void disable_ALE()
Final override for disable ALE.
unsigned u_index_axi_adv_diff() const
Overload the index at which the temperature and solute concentration variables are stored.
void get_wind_axi_adv_diff(const unsigned &ipt, const Vector< double > &s, const Vector< double > &x, Vector< double > &wind) const
Overload the wind function in the advection-diffusion equations. This provides the coupling from the ...
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 fill_in_contribution_to_jacobian(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Compute the element's residual vector and the Jacobian matrix. Jacobian is computed by finite-differe...
void fill_in_off_diagonal_jacobian_blocks_by_fd(Vector< double > &residuals, DenseMatrix< double > &jacobian)
Helper function to get the off-diagonal blocks of the Jacobian matrix by finite differences.
void enable_ALE()
Final override for enable ALE.
void compute_norm(double &norm)
std::string scalar_name_paraview(const unsigned &i) const
void output(FILE *file_pt)
C-style output function: Broken default.
AxisymmetricQAdvectionCrouzeixRaviartElement()
Constructor: call the underlying constructors and initialise the pointer to the Rayleigh number to po...
void output(FILE *file_pt, const unsigned &n_plot)
C-style output function: Broken default.
unsigned nscalar_paraview() const
Number of scalars/fields output by this element. Reimplements broken virtual function in base class.