InterfaceProblem< ELEMENT, TIMESTEPPER > Class Template Reference

Single axisymmetric fluid interface problem in rectangular domain. More...

+ Inheritance diagram for InterfaceProblem< ELEMENT, TIMESTEPPER >:

Public Member Functions

 InterfaceProblem (const unsigned &n_r, const unsigned &n_z, const double &l_r, const double &h)
 Constructor for single fluid interface problem. More...
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void actions_before_newton_convergence_check ()
 
void actions_before_newton_solve ()
 
void actions_after_newton_solve ()
 
void set_initial_condition ()
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
 InterfaceProblem (const unsigned &n_r, const unsigned &n_z, const double &l_r, const double &h)
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void actions_before_newton_convergence_check ()
 
void actions_before_newton_solve ()
 
void actions_after_newton_solve ()
 
void set_initial_condition ()
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
 InterfaceProblem ()
 Constructor. More...
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void set_initial_condition ()
 Set initial conditions. More...
 
void set_boundary_conditions ()
 Set boundary conditions. More...
 
void doc_solution (DocInfo &doc_info)
 Document the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
 InterfaceProblem (const unsigned &n_r, const unsigned &n_z1, const unsigned &n_z2, const double &l_r, const double &h1, const double &h2)
 Constructor for axisymmetric two fluid interface problem. More...
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void set_initial_condition ()
 Set initial conditions. More...
 
void set_boundary_conditions ()
 Set boundary conditions. More...
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
 InterfaceProblem (const unsigned &n_r, const unsigned &n_y, const unsigned &n_theta, const double &r_min, const double &r_max, const double &l_y, const double &theta_max)
 Problem constructor. More...
 
void actions_before_newton_convergence_check ()
 
void unsteady_run (const unsigned &nstep)
 Run an unsteady simulation with specified number of steps. More...
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
double compute_total_mass ()
 Compute the total mass. More...
 
 InterfaceProblem (const unsigned &n_r, const unsigned &n_z, const double &l_z)
 Constructor for single fluid interface problem. More...
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void actions_before_newton_convergence_check ()
 
void set_initial_condition ()
 
double global_temporal_error_norm ()
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
double compute_total_mass ()
 Compute the total mass of the insoluble surfactant. More...
 
 InterfaceProblem (const unsigned &n_r, const unsigned &n_z, const double &l_z)
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void set_initial_condition ()
 
double global_temporal_error_norm ()
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
double compute_total_mass ()
 Compute the total mass of the insoluble surfactant. More...
 
 InterfaceProblem (const unsigned &n_x, const unsigned &n_y, const double &l_x, const double &h)
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void set_initial_condition ()
 Set initial conditions. More...
 
void set_boundary_conditions ()
 Set boundary conditions. More...
 
void doc_solution (DocInfo &doc_info)
 Document the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
 InterfaceProblem (const unsigned &n_x, const unsigned &n_y, const double &l_x, const double &h)
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void set_initial_condition ()
 Set initial conditions. More...
 
void set_boundary_conditions ()
 Set boundary conditions. More...
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
 InterfaceProblem (const unsigned &Nx, const unsigned &Ny, const unsigned &Nz, const double &Lx, const double &Ly, const double &h)
 
void actions_before_newton_convergence_check ()
 
void unsteady_run (const unsigned &nstep)
 Run an unsteady simulation with specified number of steps. More...
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
 InterfaceProblem ()
 Constructor. More...
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void set_initial_condition ()
 Set initial conditions. More...
 
void set_boundary_conditions ()
 Set boundary conditions. More...
 
void doc_solution (DocInfo &doc_info)
 Document the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
 InterfaceProblem (const unsigned &n_x, const unsigned &n_y1, const unsigned &n_y2, const double &l_x, const double &h1, const double &h2)
 
 ~InterfaceProblem ()
 Destructor (empty) More...
 
void set_initial_condition ()
 Set initial conditions. More...
 
void set_boundary_conditions ()
 Set boundary conditions. More...
 
void doc_solution (DocInfo &doc_info)
 Doc the solution. More...
 
void unsteady_run (const double &t_max, const double &dt)
 Do unsteady run up to maximum time t_max with given timestep dt. More...
 
- Public Member Functions inherited from oomph::Problem
virtual void debug_hook_fct (const unsigned &i)
 
void set_analytic_dparameter (double *const &parameter_pt)
 
void unset_analytic_dparameter (double *const &parameter_pt)
 
bool is_dparameter_calculated_analytically (double *const &parameter_pt)
 
void set_analytic_hessian_products ()
 
void unset_analytic_hessian_products ()
 
bool are_hessian_products_calculated_analytically ()
 
void set_pinned_values_to_zero ()
 
bool distributed () const
 
OomphCommunicatorcommunicator_pt ()
 access function to the oomph-lib communicator More...
 
const OomphCommunicatorcommunicator_pt () const
 access function to the oomph-lib communicator, const version More...
 
 Problem ()
 
 Problem (const Problem &dummy)=delete
 Broken copy constructor. More...
 
void operator= (const Problem &)=delete
 Broken assignment operator. More...
 
virtual ~Problem ()
 Virtual destructor to clean up memory. More...
 
Mesh *& mesh_pt ()
 Return a pointer to the global mesh. More...
 
Mesh *const & mesh_pt () const
 Return a pointer to the global mesh (const version) More...
 
Mesh *& mesh_pt (const unsigned &imesh)
 
Mesh *const & mesh_pt (const unsigned &imesh) const
 Return a pointer to the i-th submesh (const version) More...
 
unsigned nsub_mesh () const
 Return number of submeshes. More...
 
unsigned add_sub_mesh (Mesh *const &mesh_pt)
 
void flush_sub_meshes ()
 
void build_global_mesh ()
 
void rebuild_global_mesh ()
 
LinearSolver *& linear_solver_pt ()
 Return a pointer to the linear solver object. More...
 
LinearSolver *const & linear_solver_pt () const
 Return a pointer to the linear solver object (const version) More...
 
LinearSolver *& mass_matrix_solver_for_explicit_timestepper_pt ()
 
LinearSolvermass_matrix_solver_for_explicit_timestepper_pt () const
 
EigenSolver *& eigen_solver_pt ()
 Return a pointer to the eigen solver object. More...
 
EigenSolver *const & eigen_solver_pt () const
 Return a pointer to the eigen solver object (const version) More...
 
Time *& time_pt ()
 Return a pointer to the global time object. More...
 
Timetime_pt () const
 Return a pointer to the global time object (const version). More...
 
doubletime ()
 Return the current value of continuous time. More...
 
double time () const
 Return the current value of continuous time (const version) More...
 
TimeStepper *& time_stepper_pt ()
 
const TimeSteppertime_stepper_pt () const
 
TimeStepper *& time_stepper_pt (const unsigned &i)
 Return a pointer to the i-th timestepper. More...
 
ExplicitTimeStepper *& explicit_time_stepper_pt ()
 Return a pointer to the explicit timestepper. More...
 
unsigned long set_timestepper_for_all_data (TimeStepper *const &time_stepper_pt, const bool &preserve_existing_data=false)
 
virtual void shift_time_values ()
 Shift all values along to prepare for next timestep. More...
 
AssemblyHandler *& assembly_handler_pt ()
 Return a pointer to the assembly handler object. More...
 
AssemblyHandler *const & assembly_handler_pt () const
 Return a pointer to the assembly handler object (const version) More...
 
doubleminimum_dt ()
 Access function to min timestep in adaptive timestepping. More...
 
doublemaximum_dt ()
 Access function to max timestep in adaptive timestepping. More...
 
unsignedmax_newton_iterations ()
 Access function to max Newton iterations before giving up. More...
 
void problem_is_nonlinear (const bool &prob_lin)
 Access function to Problem_is_nonlinear. More...
 
doublemax_residuals ()
 
booltime_adaptive_newton_crash_on_solve_fail ()
 Access function for Time_adaptive_newton_crash_on_solve_fail. More...
 
doublenewton_solver_tolerance ()
 
void add_time_stepper_pt (TimeStepper *const &time_stepper_pt)
 
void set_explicit_time_stepper_pt (ExplicitTimeStepper *const &explicit_time_stepper_pt)
 
void initialise_dt (const double &dt)
 
void initialise_dt (const Vector< double > &dt)
 
Data *& global_data_pt (const unsigned &i)
 Return a pointer to the the i-th global data object. More...
 
void add_global_data (Data *const &global_data_pt)
 
void flush_global_data ()
 
LinearAlgebraDistribution *const & dof_distribution_pt () const
 Return the pointer to the dof distribution (read-only) More...
 
unsigned long ndof () const
 Return the number of dofs. More...
 
unsigned ntime_stepper () const
 Return the number of time steppers. More...
 
unsigned nglobal_data () const
 Return the number of global data values. More...
 
unsigned self_test ()
 Self-test: Check meshes and global data. Return 0 for OK. More...
 
void enable_store_local_dof_pt_in_elements ()
 
void disable_store_local_dof_pt_in_elements ()
 
unsigned long assign_eqn_numbers (const bool &assign_local_eqn_numbers=true)
 
void describe_dofs (std::ostream &out= *(oomph_info.stream_pt())) const
 
void enable_discontinuous_formulation ()
 
void disable_discontinuous_formulation ()
 
void get_dofs (DoubleVector &dofs) const
 
void get_dofs (const unsigned &t, DoubleVector &dofs) const
 Return vector of the t'th history value of all dofs. More...
 
void set_dofs (const DoubleVector &dofs)
 Set the values of the dofs. More...
 
void set_dofs (const unsigned &t, DoubleVector &dofs)
 Set the history values of the dofs. More...
 
void set_dofs (const unsigned &t, Vector< double * > &dof_pt)
 
void add_to_dofs (const double &lambda, const DoubleVector &increment_dofs)
 Add lambda x incremenet_dofs[l] to the l-th dof. More...
 
doubleglobal_dof_pt (const unsigned &i)
 
doubledof (const unsigned &i)
 i-th dof in the problem More...
 
double dof (const unsigned &i) const
 i-th dof in the problem (const version) More...
 
double *& dof_pt (const unsigned &i)
 Pointer to i-th dof in the problem. More...
 
doubledof_pt (const unsigned &i) const
 Pointer to i-th dof in the problem (const version) More...
 
virtual void get_inverse_mass_matrix_times_residuals (DoubleVector &Mres)
 
virtual void get_dvaluesdt (DoubleVector &f)
 
virtual void get_residuals (DoubleVector &residuals)
 Get the total residuals Vector for the problem. More...
 
virtual void get_jacobian (DoubleVector &residuals, DenseDoubleMatrix &jacobian)
 
virtual void get_jacobian (DoubleVector &residuals, CRDoubleMatrix &jacobian)
 
virtual void get_jacobian (DoubleVector &residuals, CCDoubleMatrix &jacobian)
 
virtual void get_jacobian (DoubleVector &residuals, SumOfMatrices &jacobian)
 
void get_fd_jacobian (DoubleVector &residuals, DenseMatrix< double > &jacobian)
 Get the full Jacobian by finite differencing. More...
 
void get_derivative_wrt_global_parameter (double *const &parameter_pt, DoubleVector &result)
 
void get_hessian_vector_products (DoubleVectorWithHaloEntries const &Y, Vector< DoubleVectorWithHaloEntries > const &C, Vector< DoubleVectorWithHaloEntries > &product)
 
void solve_eigenproblem (const unsigned &n_eval, Vector< std::complex< double >> &eigenvalue, Vector< DoubleVector > &eigenvector, const bool &steady=true)
 Solve the eigenproblem. More...
 
void solve_eigenproblem (const unsigned &n_eval, Vector< std::complex< double >> &eigenvalue, const bool &steady=true)
 
virtual void get_eigenproblem_matrices (CRDoubleMatrix &mass_matrix, CRDoubleMatrix &main_matrix, const double &shift=0.0)
 
void assign_eigenvector_to_dofs (DoubleVector &eigenvector)
 Assign the eigenvector passed to the function to the dofs. More...
 
void add_eigenvector_to_dofs (const double &epsilon, const DoubleVector &eigenvector)
 
void store_current_dof_values ()
 Store the current values of the degrees of freedom. More...
 
void restore_dof_values ()
 Restore the stored values of the degrees of freedom. More...
 
void enable_jacobian_reuse ()
 
void disable_jacobian_reuse ()
 Disable recycling of Jacobian in Newton iteration. More...
 
bool jacobian_reuse_is_enabled ()
 Is recycling of Jacobian in Newton iteration enabled? More...
 
booluse_predictor_values_as_initial_guess ()
 
void newton_solve ()
 Use Newton method to solve the problem. More...
 
void enable_globally_convergent_newton_method ()
 enable globally convergent Newton method More...
 
void disable_globally_convergent_newton_method ()
 disable globally convergent Newton method More...
 
void newton_solve (unsigned const &max_adapt)
 
void steady_newton_solve (unsigned const &max_adapt=0)
 
void copy (Problem *orig_problem_pt)
 
virtual Problemmake_copy ()
 
virtual void read (std::ifstream &restart_file, bool &unsteady_restart)
 
virtual void read (std::ifstream &restart_file)
 
virtual void dump (std::ofstream &dump_file) const
 
void dump (const std::string &dump_file_name) const
 
void delete_all_external_storage ()
 
virtual void symmetrise_eigenfunction_for_adaptive_pitchfork_tracking ()
 
doublebifurcation_parameter_pt () const
 
void get_bifurcation_eigenfunction (Vector< DoubleVector > &eigenfunction)
 
void activate_fold_tracking (double *const &parameter_pt, const bool &block_solve=true)
 
void activate_bifurcation_tracking (double *const &parameter_pt, const DoubleVector &eigenvector, const bool &block_solve=true)
 
void activate_bifurcation_tracking (double *const &parameter_pt, const DoubleVector &eigenvector, const DoubleVector &normalisation, const bool &block_solve=true)
 
void activate_pitchfork_tracking (double *const &parameter_pt, const DoubleVector &symmetry_vector, const bool &block_solve=true)
 
void activate_hopf_tracking (double *const &parameter_pt, const bool &block_solve=true)
 
void activate_hopf_tracking (double *const &parameter_pt, const double &omega, const DoubleVector &null_real, const DoubleVector &null_imag, const bool &block_solve=true)
 
void deactivate_bifurcation_tracking ()
 
void reset_assembly_handler_to_default ()
 Reset the system to the standard non-augemented state. More...
 
double arc_length_step_solve (double *const &parameter_pt, const double &ds, const unsigned &max_adapt=0)
 
double arc_length_step_solve (Data *const &data_pt, const unsigned &data_index, const double &ds, const unsigned &max_adapt=0)
 
void reset_arc_length_parameters ()
 
intsign_of_jacobian ()
 
void explicit_timestep (const double &dt, const bool &shift_values=true)
 Take an explicit timestep of size dt. More...
 
void unsteady_newton_solve (const double &dt)
 
void unsteady_newton_solve (const double &dt, const bool &shift_values)
 
void unsteady_newton_solve (const double &dt, const unsigned &max_adapt, const bool &first, const bool &shift=true)
 
double doubly_adaptive_unsteady_newton_solve (const double &dt, const double &epsilon, const unsigned &max_adapt, const bool &first, const bool &shift=true)
 
double doubly_adaptive_unsteady_newton_solve (const double &dt, const double &epsilon, const unsigned &max_adapt, const unsigned &suppress_resolve_after_spatial_adapt_flag, const bool &first, const bool &shift=true)
 
double adaptive_unsteady_newton_solve (const double &dt_desired, const double &epsilon)
 
double adaptive_unsteady_newton_solve (const double &dt_desired, const double &epsilon, const bool &shift_values)
 
void assign_initial_values_impulsive ()
 
void assign_initial_values_impulsive (const double &dt)
 
void calculate_predictions ()
 Calculate predictions. More...
 
void enable_mass_matrix_reuse ()
 
void disable_mass_matrix_reuse ()
 
bool mass_matrix_reuse_is_enabled ()
 Return whether the mass matrix is being reused. More...
 
void refine_uniformly (const Vector< unsigned > &nrefine_for_mesh)
 
void refine_uniformly (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh)
 
void refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void refine_uniformly (DocInfo &doc_info)
 
void refine_uniformly_and_prune (DocInfo &doc_info)
 
void refine_uniformly ()
 
void refine_uniformly (const unsigned &i_mesh, DocInfo &doc_info)
 Do uniform refinement for submesh i_mesh with documentation. More...
 
void refine_uniformly (const unsigned &i_mesh)
 Do uniform refinement for submesh i_mesh without documentation. More...
 
void p_refine_uniformly (const Vector< unsigned > &nrefine_for_mesh)
 
void p_refine_uniformly (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void p_refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh)
 
void p_refine_uniformly_and_prune (const Vector< unsigned > &nrefine_for_mesh, DocInfo &doc_info)
 
void p_refine_uniformly (DocInfo &doc_info)
 
void p_refine_uniformly_and_prune (DocInfo &doc_info)
 
void p_refine_uniformly ()
 
void p_refine_uniformly (const unsigned &i_mesh, DocInfo &doc_info)
 Do uniform p-refinement for submesh i_mesh with documentation. More...
 
void p_refine_uniformly (const unsigned &i_mesh)
 Do uniform p-refinement for submesh i_mesh without documentation. More...
 
void refine_selected_elements (const Vector< unsigned > &elements_to_be_refined)
 
void refine_selected_elements (const Vector< RefineableElement * > &elements_to_be_refined_pt)
 
void refine_selected_elements (const unsigned &i_mesh, const Vector< unsigned > &elements_to_be_refined)
 
void refine_selected_elements (const unsigned &i_mesh, const Vector< RefineableElement * > &elements_to_be_refined_pt)
 
void refine_selected_elements (const Vector< Vector< unsigned >> &elements_to_be_refined)
 
void refine_selected_elements (const Vector< Vector< RefineableElement * >> &elements_to_be_refined_pt)
 
void p_refine_selected_elements (const Vector< unsigned > &elements_to_be_refined)
 
void p_refine_selected_elements (const Vector< PRefineableElement * > &elements_to_be_refined_pt)
 
void p_refine_selected_elements (const unsigned &i_mesh, const Vector< unsigned > &elements_to_be_refined)
 
void p_refine_selected_elements (const unsigned &i_mesh, const Vector< PRefineableElement * > &elements_to_be_refined_pt)
 
void p_refine_selected_elements (const Vector< Vector< unsigned >> &elements_to_be_refined)
 
void p_refine_selected_elements (const Vector< Vector< PRefineableElement * >> &elements_to_be_refined_pt)
 
unsigned unrefine_uniformly ()
 
unsigned unrefine_uniformly (const unsigned &i_mesh)
 
void p_unrefine_uniformly (DocInfo &doc_info)
 
void p_unrefine_uniformly (const unsigned &i_mesh, DocInfo &doc_info)
 Do uniform p-unrefinement for submesh i_mesh without documentation. More...
 
void adapt (unsigned &n_refined, unsigned &n_unrefined)
 
void adapt ()
 
void p_adapt (unsigned &n_refined, unsigned &n_unrefined)
 
void p_adapt ()
 
void adapt_based_on_error_estimates (unsigned &n_refined, unsigned &n_unrefined, Vector< Vector< double >> &elemental_error)
 
void adapt_based_on_error_estimates (Vector< Vector< double >> &elemental_error)
 
void get_all_error_estimates (Vector< Vector< double >> &elemental_error)
 
void doc_errors (DocInfo &doc_info)
 Get max and min error for all elements in submeshes. More...
 
void doc_errors ()
 Get max and min error for all elements in submeshes. More...
 
void enable_info_in_newton_solve ()
 
void disable_info_in_newton_solve ()
 Disable the output of information when in the newton solver. More...
 
- Public Member Functions inherited from oomph::ExplicitTimeSteppableObject
 ExplicitTimeSteppableObject ()
 Empty constructor. More...
 
 ExplicitTimeSteppableObject (const ExplicitTimeSteppableObject &)=delete
 Broken copy constructor. More...
 
void operator= (const ExplicitTimeSteppableObject &)=delete
 Broken assignment operator. More...
 
virtual ~ExplicitTimeSteppableObject ()
 Empty destructor. More...
 
virtual void actions_before_explicit_stage ()
 
virtual void actions_after_explicit_stage ()
 

Public Attributes

HorizontalSingleLayerSpineMesh< ELEMENT > * Bulk_mesh_pt
 Access function for the specific mesh. More...
 
MeshInterface_mesh_pt
 Mesh for the free surface (interface) elements. More...
 
SingleLayerSpineMesh< ELEMENT > * Bulk_mesh_pt
 Access function for the specific mesh. More...
 
MeshSurface_mesh_pt
 Mesh for the interface elements. More...
 
MyHorizontalSingleLayerSpineMesh< ELEMENT > * Bulk_mesh_pt
 Access function for the specific mesh. More...
 
ElasticRectangularQuadMesh< ELEMENT > * Bulk_mesh_pt
 Access function for the specific mesh. More...
 
- Public Attributes inherited from oomph::Problem
bool Shut_up_in_newton_solve
 

Private Member Functions

void deform_free_surface (const double &epsilon)
 Deform the mesh/free surface to a prescribed function. More...
 
void deform_free_surface (const double &epsilon, const double &k)
 Deform the mesh/free surface to a prescribed function. More...
 
void actions_before_newton_solve ()
 No actions required before solve step. More...
 
void actions_after_newton_solve ()
 No actions required after solve step. More...
 
void actions_before_implicit_timestep ()
 
void actions_before_adapt ()
 Strip off the interface elements before adapting the bulk mesh. More...
 
void actions_after_adapt ()
 Rebuild the mesh of interface elements after adapting the bulk mesh. More...
 
void create_interface_elements ()
 Create the 1d interface elements. More...
 
void delete_interface_elements ()
 Delete the 1d interface elements. More...
 
void deform_free_surface (const double &epsilon, const double &k)
 Deform the mesh/free surface to a prescribed function. More...
 
void fix_pressure (const unsigned &e, const unsigned &pdof, const double &pvalue)
 Fix pressure in element e at pressure dof pdof and set to pvalue. More...
 
void actions_before_newton_convergence_check ()
 
void actions_before_newton_solve ()
 No actions required before solve step. More...
 
void actions_after_newton_solve ()
 
void deform_free_surface (const double &epsilon, const double &k)
 Deform the mesh/free surface to a prescribed function. More...
 
void fix_pressure (const unsigned &e, const unsigned &pdof, const double &pvalue)
 Fix pressure in element e at pressure dof pdof and set to pvalue. More...
 
void deform_free_surface (const double &epsilon)
 Deform the mesh/free surface to a prescribed function. More...
 
void deform_free_surface (const double &epsilon)
 Deform the mesh/free surface to a prescribed function. More...
 
void actions_before_newton_solve ()
 No actions required before solve step. More...
 
void actions_after_newton_solve ()
 No actions required after solve step. More...
 
void actions_before_implicit_timestep ()
 
void deform_free_surface (const double &epsilon, const unsigned &n_periods)
 Deform the mesh/free surface to a prescribed function. More...
 
void actions_before_newton_convergence_check ()
 
void actions_before_newton_solve ()
 No actions required before solve step. More...
 
void actions_after_newton_solve ()
 
void deform_free_surface (const double &epsilon, const unsigned &n_periods)
 Deform the mesh/free surface to a prescribed function. More...
 
void actions_before_newton_solve ()
 No actions required before solve step. More...
 
void actions_after_newton_solve ()
 No actions required after solve step. More...
 
void actions_before_implicit_timestep ()
 
void actions_before_adapt ()
 Strip off the interface elements before adapting the bulk mesh. More...
 
void actions_after_adapt ()
 Rebuild the mesh of interface elements after adapting the bulk mesh. More...
 
void create_interface_elements ()
 Create the 1d interface elements. More...
 
void delete_interface_elements ()
 Delete the 1d interface elements. More...
 
void deform_free_surface (const double &epsilon, const unsigned &n_periods)
 Deform the mesh/free surface to a prescribed function. More...
 
void fix_pressure (const unsigned &e, const unsigned &pdof, const double &pvalue)
 Fix pressure in element e at pressure dof pdof and set to pvalue. More...
 
void actions_before_newton_convergence_check ()
 
void actions_before_newton_solve ()
 No actions required before solve step. More...
 
void actions_after_newton_solve ()
 
void deform_free_surface (const double &epsilon, const unsigned &n_periods)
 Deform the mesh/free surface to a prescribed function. More...
 
void fix_pressure (const unsigned &e, const unsigned &pdof, const double &pvalue)
 Fix pressure in element e at pressure dof pdof and set to pvalue. More...
 

Private Attributes

ofstream Trace_file
 Trace file. More...
 
double Lr
 Width of domain. More...
 
double Height
 Height of the domain. More...
 
ElasticRefineableTwoLayerMesh< ELEMENT > * Bulk_mesh_pt
 Pointer to the (specific) "bulk" mesh. More...
 
ConstitutiveLawConstitutive_law_pt
 Pointer to the constitutive law. More...
 
TwoLayerSpineMesh< ELEMENT > * Bulk_mesh_pt
 Pointer to the specific bulk mesh. More...
 
double R_max
 Axial lengths of domain. More...
 
double L_y
 
AnnularSpineMesh< ELEMENT > * Bulk_mesh_pt
 Pointer to bulk mesh. More...
 
NodeDocument_node_pt
 Pointer to a node for documentation purposes. More...
 
double Lx
 Width of domain. More...
 
double Ly
 
SingleLayerCubicSpineMesh< ELEMENT > * Bulk_mesh_pt
 Pointer to bulk mesh. More...
 

Additional Inherited Members

- Public Types inherited from oomph::Problem
typedef void(* SpatialErrorEstimatorFctPt) (Mesh *&mesh_pt, Vector< double > &elemental_error)
 Function pointer for spatial error estimator. More...
 
typedef void(* SpatialErrorEstimatorWithDocFctPt) (Mesh *&mesh_pt, Vector< double > &elemental_error, DocInfo &doc_info)
 Function pointer for spatial error estimator with doc. More...
 
- Static Public Attributes inherited from oomph::Problem
static bool Suppress_warning_about_actions_before_read_unstructured_meshes
 
- Protected Types inherited from oomph::Problem
enum  Assembly_method {
  Perform_assembly_using_vectors_of_pairs , Perform_assembly_using_two_vectors , Perform_assembly_using_maps , Perform_assembly_using_lists ,
  Perform_assembly_using_two_arrays
}
 Enumerated flags to determine which sparse assembly method is used. More...
 
- Protected Member Functions inherited from oomph::Problem
unsigned setup_element_count_per_dof ()
 
virtual void sparse_assemble_row_or_column_compressed (Vector< int * > &column_or_row_index, Vector< int * > &row_or_column_start, Vector< double * > &value, Vector< unsigned > &nnz, Vector< double * > &residual, bool compressed_row_flag)
 
virtual void actions_before_newton_step ()
 
virtual void actions_after_newton_step ()
 
virtual void actions_after_implicit_timestep ()
 
virtual void actions_after_implicit_timestep_and_error_estimation ()
 
virtual void actions_before_explicit_timestep ()
 Actions that should be performed before each explicit time step. More...
 
virtual void actions_after_explicit_timestep ()
 Actions that should be performed after each explicit time step. More...
 
virtual void actions_before_read_unstructured_meshes ()
 
virtual void actions_after_read_unstructured_meshes ()
 
virtual void actions_after_change_in_global_parameter (double *const &parameter_pt)
 
virtual void actions_after_change_in_bifurcation_parameter ()
 
virtual void actions_after_parameter_increase (double *const &parameter_pt)
 
doubledof_derivative (const unsigned &i)
 
doubledof_current (const unsigned &i)
 
unsigned newton_solve_continuation (double *const &parameter_pt)
 
unsigned newton_solve_continuation (double *const &parameter_pt, DoubleVector &z)
 
void calculate_continuation_derivatives (double *const &parameter_pt)
 
void calculate_continuation_derivatives (const DoubleVector &z)
 
void calculate_continuation_derivatives_fd (double *const &parameter_pt)
 
bool does_pointer_correspond_to_problem_data (double *const &parameter_pt)
 
void set_consistent_pinned_values_for_continuation ()
 
- Protected Attributes inherited from oomph::Problem
Vector< Problem * > Copy_of_problem_pt
 
std::map< double *, boolCalculate_dparameter_analytic
 
bool Calculate_hessian_products_analytic
 
LinearAlgebraDistributionDof_distribution_pt
 
Vector< double * > Dof_pt
 Vector of pointers to dofs. More...
 
DoubleVectorWithHaloEntries Element_count_per_dof
 
double Relaxation_factor
 
double Newton_solver_tolerance
 
unsigned Max_newton_iterations
 Maximum number of Newton iterations. More...
 
unsigned Nnewton_iter_taken
 
Vector< doubleMax_res
 Maximum residuals at start and after each newton iteration. More...
 
double Max_residuals
 
bool Time_adaptive_newton_crash_on_solve_fail
 
bool Jacobian_reuse_is_enabled
 Is re-use of Jacobian in Newton iteration enabled? Default: false. More...
 
bool Jacobian_has_been_computed
 
bool Problem_is_nonlinear
 
bool Pause_at_end_of_sparse_assembly
 
bool Doc_time_in_distribute
 
unsigned Sparse_assembly_method
 
unsigned Sparse_assemble_with_arrays_initial_allocation
 
unsigned Sparse_assemble_with_arrays_allocation_increment
 
Vector< Vector< unsigned > > Sparse_assemble_with_arrays_previous_allocation
 
double Numerical_zero_for_sparse_assembly
 
double FD_step_used_in_get_hessian_vector_products
 
bool Mass_matrix_reuse_is_enabled
 
bool Mass_matrix_has_been_computed
 
bool Discontinuous_element_formulation
 
double Minimum_dt
 Minimum desired dt: if dt falls below this value, exit. More...
 
double Maximum_dt
 Maximum desired dt. More...
 
double DTSF_max_increase
 
double DTSF_min_decrease
 
double Minimum_dt_but_still_proceed
 
bool Scale_arc_length
 Boolean to control whether arc-length should be scaled. More...
 
double Desired_proportion_of_arc_length
 Proportion of the arc-length to taken by the parameter. More...
 
double Theta_squared
 
int Sign_of_jacobian
 Storage for the sign of the global Jacobian. More...
 
double Continuation_direction
 
double Parameter_derivative
 Storage for the derivative of the global parameter wrt arc-length. More...
 
double Parameter_current
 Storage for the present value of the global parameter. More...
 
bool Use_continuation_timestepper
 Boolean to control original or new storage of dof stuff. More...
 
unsigned Dof_derivative_offset
 
unsigned Dof_current_offset
 
Vector< doubleDof_derivative
 Storage for the derivative of the problem variables wrt arc-length. More...
 
Vector< doubleDof_current
 Storage for the present values of the variables. More...
 
double Ds_current
 Storage for the current step value. More...
 
unsigned Desired_newton_iterations_ds
 
double Minimum_ds
 Minimum desired value of arc-length. More...
 
bool Bifurcation_detection
 Boolean to control bifurcation detection via determinant of Jacobian. More...
 
bool Bisect_to_find_bifurcation
 Boolean to control wheter bisection is used to located bifurcation. More...
 
bool First_jacobian_sign_change
 Boolean to indicate whether a sign change has occured in the Jacobian. More...
 
bool Arc_length_step_taken
 Boolean to indicate whether an arc-length step has been taken. More...
 
bool Use_finite_differences_for_continuation_derivatives
 
OomphCommunicatorCommunicator_pt
 The communicator for this problem. More...
 
bool Always_take_one_newton_step
 
double Timestep_reduction_factor_after_nonconvergence
 
bool Keep_temporal_error_below_tolerance
 
- Static Protected Attributes inherited from oomph::Problem
static ContinuationStorageScheme Continuation_time_stepper
 Storage for the single static continuation timestorage object. More...
 

Detailed Description

template<class ELEMENT, class TIMESTEPPER>
class InterfaceProblem< ELEMENT, TIMESTEPPER >

Single axisymmetric fluid interface problem in rectangular domain.

Single fluid interface problem.

Axisymmetric two fluid interface problem in a rectangular domain.

/////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////////

Single fluid interface problem including transport of an insoluble surfactant.

Single axisymmetric fluid interface problem including the transport of an insoluble surfactant.

Single axisymmetric fluid interface problem including the transport of an soluble surfactant.

/////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////// Single fluid free surface problem in a rectangular domain which is periodic in the x direction

/////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////// Two fluid interface problem in a rectangular domain which is periodic in the x direction

Constructor & Destructor Documentation

◆ InterfaceProblem() [1/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_r,
const unsigned n_z,
const double l_r,
const double h 
)

Constructor for single fluid interface problem.

Constructor for single fluid free surface problem.

Constructor: Pass the number of elements and the lengths of the domain in the r and z directions (h is the height of the fluid layer i.e. the length of the domain in the z direction)

203  : Lr(l_r),
204  Height(h)
205 
206 {
207 
208  // Allocate the timestepper (this constructs the time object as well)
209  add_time_stepper_pt(new TIMESTEPPER);
210 
211  // Build and assign mesh (the "false" boolean flag tells the mesh
212  // constructor that the domain is not periodic in r)
213  Bulk_mesh_pt =
215 
216  //Create "surface mesh" that will only contain the interface elements
217  Interface_mesh_pt = new Mesh;
218  {
219  // How many bulk elements are adjacent to boundary b?
220  unsigned n_element = Bulk_mesh_pt->nboundary_element(1);
221 
222  // Loop over the bulk elements adjacent to boundary b?
223  for(unsigned e=0;e<n_element;e++)
224  {
225  // Get pointer to the bulk element that is adjacent to boundary b
226  ELEMENT* bulk_elem_pt = dynamic_cast<ELEMENT*>(
227  Bulk_mesh_pt->boundary_element_pt(1,e));
228 
229  // Find the index of the face of element e along boundary b
230  int face_index = Bulk_mesh_pt->face_index_at_boundary(1,e);
231 
232  // Build the corresponding free surface element
233  SpineAxisymmetricFluidInterfaceElement<ELEMENT>* interface_element_pt = new
234  SpineAxisymmetricFluidInterfaceElement<ELEMENT>(bulk_elem_pt,face_index);
235 
236  //Add the prescribed-flux element to the surface mesh
237  Interface_mesh_pt->add_element_pt(interface_element_pt);
238 
239  } //end of loop over bulk elements adjacent to boundary b
240  }
241 
242 
243  // Add the two sub meshes to the problem
246 
247  // Combine all submeshes into a single Mesh
249 
250 
251  // --------------------------------------------
252  // Set the boundary conditions for this problem
253  // --------------------------------------------
254 
255  // All nodes are free by default -- just pin the ones that have
256  // Dirichlet conditions here
257 
258  // Determine number of mesh boundaries
259  const unsigned n_boundary = Bulk_mesh_pt->nboundary();
260 
261  // Loop over mesh boundaries
262  for(unsigned b=0;b<n_boundary;b++)
263  {
264  // Determine number of nodes on boundary b
265  const unsigned n_node = Bulk_mesh_pt->nboundary_node(b);
266 
267  // Loop over nodes on boundary b
268  for(unsigned n=0;n<n_node;n++)
269  {
270  // Pin azimuthal velocity on bounds
271  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(2);
272 
273  // Pin radial velocity on axis of symmetry
274  if(b==3)
275  {
276  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(0);
277  }
278  // Pin axial velocity on bottom and top walls (no penetration)
279  if(b==0 || b==2)
280  {
281  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(1);
282  }
283  } // End of loop over nodes on boundary b
284  } // End of loop over mesh boundaries
285 
286  // ----------------------------------------------------------------
287  // Complete the problem setup to make the elements fully functional
288  // ----------------------------------------------------------------
289 
290  // Determine number of bulk elements in mesh
291  const unsigned n_bulk = Bulk_mesh_pt->nelement();
292 
293  // Loop over the bulk elements
294  for(unsigned e=0;e<n_bulk;e++)
295  {
296  // Upcast from GeneralisedElement to the present element
297  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->element_pt(e));
298 
299  // Set the Reynolds number
300  el_pt->re_pt() = &Global_Physical_Variables::Re;
301 
302  // Set the Womersley number
303  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
304 
305  // Set the product of the Reynolds number and the inverse of the
306  // Froude number
307  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
308 
309  // Set the direction of gravity
310  el_pt->g_pt() = &Global_Physical_Variables::G;
311 
312  } // End of loop over bulk elements
313 
314  // Create a Data object whose single value stores the external pressure
315  Data* external_pressure_data_pt = new Data(1);
316 
317  // Pin and set the external pressure to some arbitrary value
318  double p_ext = Global_Physical_Variables::P_ext;
319 
320  external_pressure_data_pt->pin(0);
321  external_pressure_data_pt->set_value(0,p_ext);
322 
323  // Determine number of 1D interface elements in mesh
324  const unsigned n_interface_element = Interface_mesh_pt->nelement();
325 
326  // Loop over the interface elements
327  for(unsigned e=0;e<n_interface_element;e++)
328  {
329  // Upcast from GeneralisedElement to the present element
333 
334  // Set the Capillary number
336 
337  // Pass the Data item that contains the single external pressure value
338  el_pt->set_external_pressure_data(external_pressure_data_pt);
339 
340  } // End of loop over interface elements
341 
342  // Setup equation numbering scheme
343  cout << "Number of equations: " << assign_eqn_numbers() << std::endl;
344 
345 } // End of constructor
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
Array< double, 1, 3 > e(1./3., 0.5, 2.)
Scalar * b
Definition: benchVecAdd.cpp:17
double Lr
Width of domain.
Definition: rayleigh_instability.cc:186
double Height
Height of the domain.
Definition: rayleigh_instability.cc:189
HorizontalSingleLayerSpineMesh< ELEMENT > * Bulk_mesh_pt
Access function for the specific mesh.
Definition: rayleigh_instability.cc:145
Mesh * Interface_mesh_pt
Mesh for the free surface (interface) elements.
Definition: rayleigh_instability.cc:148
Definition: nodes.h:86
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:385
void set_value(const unsigned &i, const double &value_)
Definition: nodes.h:271
double *& ca_pt()
Pointer to the Capillary number.
Definition: interface_elements.h:492
void set_external_pressure_data(Data *external_pressure_data_pt)
Definition: interface_elements.h:539
Definition: horizontal_single_layer_spine_mesh.template.h:47
Definition: mesh.h:67
GeneralisedElement *& element_pt(const unsigned long &e)
Return pointer to element e.
Definition: mesh.h:448
void add_element_pt(GeneralisedElement *const &element_pt)
Add a (pointer to) an element to the mesh.
Definition: mesh.h:617
unsigned long nelement() const
Return number of elements in the mesh.
Definition: mesh.h:590
void add_time_stepper_pt(TimeStepper *const &time_stepper_pt)
Definition: problem.cc:1545
unsigned add_sub_mesh(Mesh *const &mesh_pt)
Definition: problem.h:1330
void build_global_mesh()
Definition: problem.cc:1493
unsigned long assign_eqn_numbers(const bool &assign_local_eqn_numbers=true)
Definition: problem.cc:1989
TimeStepper *& time_stepper_pt()
Definition: problem.h:1524
Definition: specific_node_update_interface_elements.h:592
double P_ext
External pressure.
Definition: fibre.cc:64
double ReSt
Womersley number.
Definition: rayleigh_instability.cc:56
double Ca
Capillary number.
Definition: fibre.cc:61
double ReInvFr
Product of Reynolds number and inverse of Froude number.
Definition: fibre.cc:58
double Re
Reynolds number.
Definition: fibre.cc:55
Vector< double > G(3)
Direction of gravity.
Definition: spherical_shell_convection.cc:62

References oomph::Mesh::add_element_pt(), oomph::Problem::add_sub_mesh(), oomph::Problem::add_time_stepper_pt(), oomph::Problem::assign_eqn_numbers(), b, oomph::Problem::build_global_mesh(), InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt, Global_Physical_Variables::Ca, oomph::FluidInterfaceElement::ca_pt(), e(), oomph::Mesh::element_pt(), Global_Physical_Variables::G, InterfaceProblem< ELEMENT, TIMESTEPPER >::Interface_mesh_pt, n, oomph::Mesh::nelement(), Global_Physical_Variables::P_ext, oomph::Data::pin(), Global_Physical_Variables::Re, Global_Physical_Variables::ReInvFr, Global_Physical_Variables::ReSt, oomph::FluidInterfaceElement::set_external_pressure_data(), oomph::Data::set_value(), and oomph::Problem::time_stepper_pt().

◆ ~InterfaceProblem() [1/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

100 {}

◆ InterfaceProblem() [2/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_r,
const unsigned n_z,
const double l_r,
const double h 
)

Constructor: Pass the number of elements and the lengths of the domain in the r and z directions (h is the height of the fluid layer i.e. the length of the domain in the z direction)

◆ ~InterfaceProblem() [2/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

93 {}

◆ InterfaceProblem() [3/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem

Constructor.

Constructor for two fluid interface problem.

Constructor for axisymmetric two fluid interface problem.

261 {
262  // Allocate the timestepper (this constructs the time object as well)
263  add_time_stepper_pt(new TIMESTEPPER);
264 
265  // Define number of elements in r direction
266  const unsigned n_r = 3;
267 
268  // Define number of elements in z direction in lower fluid (fluid 1)
269  const unsigned n_z1 = 3;
270 
271  // Define number of elements in z direction in upper fluid (fluid 2)
272  const unsigned n_z2 = 3;
273 
274  // Define width of domain and store as class member data
275  const double l_r = 1.0;
276  this->Lr = l_r;
277 
278  // Define height of lower fluid layer
279  const double h1 = 1.0;
280 
281  // Define height of upper fluid layer
282  const double h2 = 1.0;
283 
284  // Build and assign the "bulk" mesh
286  (n_r,n_z1,n_z2,l_r,h1,h2,time_stepper_pt());
287 
288  // Create and set the error estimator for spatial adaptivity
289  Bulk_mesh_pt->spatial_error_estimator_pt() = new Z2ErrorEstimator;
290 
291  // Set the maximum refinement level for the mesh to 4
292  Bulk_mesh_pt->max_refinement_level() = 4;
293 
294  // Create the "surface" mesh that will contain only the interface
295  // elements. The constructor just creates the mesh without giving
296  // it any elements, nodes, etc.
297  Surface_mesh_pt = new Mesh;
298 
299  // Create interface elements at the boundary between the two fluids,
300  // and add them to the surface mesh
302 
303  // Add the two sub meshes to the problem
306 
307  // Combine all sub-meshes into a single mesh
309 
310  // --------------------------------------------
311  // Set the boundary conditions for this problem
312  // --------------------------------------------
313 
314  // All nodes are free by default -- just pin the ones that have
315  // Dirichlet conditions here
316 
317  // Determine number of mesh boundaries
318  const unsigned n_boundary = Bulk_mesh_pt->nboundary();
319 
320  // Loop over mesh boundaries
321  for(unsigned b=0;b<n_boundary;b++)
322  {
323  // Determine number of nodes on boundary b
324  const unsigned n_node = Bulk_mesh_pt->nboundary_node(b);
325 
326  // Loop over nodes on boundary b
327  for(unsigned n=0;n<n_node;n++)
328  {
329  // Fluid boundary conditions:
330  // --------------------------
331 
332  // Pin radial and azimuthal velocities (no slip/penetration)
333  // on all boundaries other than the interface (b=4)
334  if(b!=4)
335  {
336  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(0);
337  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(2);
338  }
339 
340  // Pin axial velocity on top (b=2) and bottom (b=0) boundaries
341  // (no penetration). Because we have a slippery outer wall we do
342  // NOT pin the axial velocity on this boundary (b=1); similarly,
343  // we do not pin the axial velocity on the symmetry boundary (b=3).
344  if(b==0 || b==2) { Bulk_mesh_pt->boundary_node_pt(b,n)->pin(1); }
345 
346  // Solid boundary conditions:
347  // --------------------------
348 
349  // Pin vertical displacement on solid boundaries
350  if(b==0 || b==2) { Bulk_mesh_pt->boundary_node_pt(b,n)->pin_position(1); }
351 
352  } // End of loop over nodes on boundary b
353  } // End of loop over mesh boundaries
354 
355  // Loop over all nodes in mesh
356  const unsigned n_node = Bulk_mesh_pt->nnode();
357  for(unsigned n=0;n<n_node;n++)
358  {
359  // Pin horizontal displacement of all nodes
360  Bulk_mesh_pt->node_pt(n)->pin_position(0);
361 
362  // Pin all azimuthal velocities throughout the bulk of the domain
363  Bulk_mesh_pt->node_pt(n)->pin(2);
364  }
365 
366  // Define a constitutive law for the solid equations: generalised Hookean
368 
369  // ----------------------------------------------------------------
370  // Complete the problem setup to make the elements fully functional
371  // ----------------------------------------------------------------
372 
373  // Compute number of bulk elements in lower/upper fluids
374  const unsigned n_lower = n_r*n_z1;
375  const unsigned n_upper = n_r*n_z2;
376 
377  // Loop over bulk elements in lower fluid
378  for(unsigned e=0;e<n_lower;e++)
379  {
380  // Upcast from GeneralisedElement to the present element
381  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->element_pt(e));
382 
383  // Set the Reynolds number
384  el_pt->re_pt() = &Global_Physical_Variables::Re;
385 
386  // Set the Womersley number
387  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
388 
389  // Set the product of the Reynolds number and the inverse of the
390  // Froude number
391  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
392 
393  // Set the direction of gravity
394  el_pt->g_pt() = &Global_Physical_Variables::G;
395 
396  // Set the constitutive law
397  el_pt->constitutive_law_pt() = Constitutive_law_pt;
398 
399  } // End of loop over bulk elements in lower fluid
400 
401  // Loop over bulk elements in upper fluid
402  for(unsigned e=n_lower;e<(n_lower+n_upper);e++)
403  {
404  // Upcast from GeneralisedElement to the present element
405  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->element_pt(e));
406 
407  // Set the Reynolds number
408  el_pt->re_pt() = &Global_Physical_Variables::Re;
409 
410  // Set the Womersley number
411  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
412 
413  // Set the product of the Reynolds number and the inverse of the
414  // Froude number
415  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
416 
417  // Set the direction of gravity
418  el_pt->g_pt() = &Global_Physical_Variables::G;
419 
420  // Set the viscosity ratio
421  el_pt->viscosity_ratio_pt() = &Global_Physical_Variables::Viscosity_Ratio;
422 
423  // Set the density ratio
424  el_pt->density_ratio_pt() = &Global_Physical_Variables::Density_Ratio;
425 
426  // Set the constitutive law
427  el_pt->constitutive_law_pt() = Constitutive_law_pt;
428 
429  } // End of loop over bulk elements in upper fluid
430 
431  // Set the pressure in the first element at 'node' 0 to 0.0
432  fix_pressure(0,0,0.0);
433 
434  // Pin the redundant solid pressures (if any)
436  Bulk_mesh_pt->element_pt());
437 
438  // Apply the boundary conditions
440 
441  // Set up equation numbering scheme
442  cout << "Number of equations: " << assign_eqn_numbers() << std::endl;
443 
444 } // End of constructor
Definition: elastic_two_layer_interface_axisym.cc:108
Mesh * Surface_mesh_pt
Mesh for the interface elements.
Definition: single_layer_free_surface_axisym.cc:140
void set_boundary_conditions()
Set boundary conditions.
Definition: elastic_two_layer_interface_axisym.cc:484
ConstitutiveLaw * Constitutive_law_pt
Pointer to the constitutive law.
Definition: elastic_two_layer_interface_axisym.cc:243
void create_interface_elements()
Create the 1d interface elements.
Definition: elastic_two_layer_interface_axisym.cc:577
void fix_pressure(const unsigned &e, const unsigned &pdof, const double &pvalue)
Fix pressure in element e at pressure dof pdof and set to pvalue.
Definition: elastic_two_layer_interface_axisym.cc:227
Definition: constitutive_laws.h:699
Definition: solid_elements.h:58
Definition: error_estimator.h:266
double Nu
Poisson's ratio.
Definition: TwenteMeshGluing.cpp:71
double Density_Ratio
Definition: elastic_two_layer_interface_axisym.cc:80
double Viscosity_Ratio
Definition: elastic_two_layer_interface_axisym.cc:76

References b, Constitutive::Constitutive_law_pt, Global_Physical_Variables::Density_Ratio, e(), Global_Physical_Variables::G, Global_Parameters::Lr, n, Global_Physical_Variables::Nu, Global_Physical_Variables::Re, Global_Physical_Variables::ReInvFr, Global_Physical_Variables::ReSt, and Global_Physical_Variables::Viscosity_Ratio.

◆ ~InterfaceProblem() [3/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

182 {}

◆ InterfaceProblem() [4/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_r,
const unsigned n_z1,
const unsigned n_z2,
const double l_r,
const double h1,
const double h2 
)

Constructor for axisymmetric two fluid interface problem.

Constructor for two fluid interface problem.

Problem constructor.

Constructor: Pass the number of elements and the width of the domain in the r direction. Also pass the number of elements in the z direction of the bottom (fluid 1) and top (fluid 2) layers, along with the heights of both layers.

178  : Lr(l_r)
179 {
180 
181  // Allocate the timestepper (this constructs the time object as well)
182  add_time_stepper_pt(new TIMESTEPPER);
183 
184  // Build and assign mesh
186  (n_r,n_z1,n_z2,l_r,h1,h2,time_stepper_pt());
187 
188  Surface_mesh_pt = new Mesh;
189 
190  //Loop over the horizontal elements
191  for(unsigned i=0;i<n_r;i++)
192  {
193  //Construct a new 1D line element on the face on which the local
194  //coordinate 1 is fixed at its max. value (1) -- Face 2
195  FiniteElement *interface_element_pt = new
197  Bulk_mesh_pt->finite_element_pt(n_r*(n_z1-1)+i),2);
198 
199  this->Surface_mesh_pt->add_element_pt(interface_element_pt);
200  }
201 
202  // --------------------------------------------
203  // Set the boundary conditions for this problem
204  // --------------------------------------------
205 
206  // All nodes are free by default -- just pin the ones that have
207  // Dirichlet conditions here
208 
209  // Loop over mesh boundaries
210  for(unsigned b=0;b<6;b++)
211  {
212  // Determine number of nodes on boundary b
213  const unsigned n_node = Bulk_mesh_pt->nboundary_node(b);
214 
215  // Loop over nodes on boundary b
216  for (unsigned n=0;n<n_node;n++)
217  {
218  // Pin radial and azimuthal velocities on all boundaries
219  // (no slip/penetration)
220  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(0);
221  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(2);
222 
223  // Pin axial velocity on top (b=3) and bottom (b=0) boundaries
224  // (no penetration). Because we have a slippery outer wall we do
225  // NOT pin the axial velocity on this boundary (b=1,2); similarly,
226  // we do not pin the axial velocity on the symmetry boundary (b=4,5).
227  if(b==0 || b==3) { Bulk_mesh_pt->boundary_node_pt(b,n)->pin(1); }
228  }
229  }
230 
231  // Pin all azimuthal velocities throughout the bulk of the domain
232  const unsigned n_node = Bulk_mesh_pt->nnode();
233  for(unsigned n=0;n<n_node;n++) { Bulk_mesh_pt->node_pt(n)->pin(2); }
234 
235  // ----------------------------------------------------------------
236  // Complete the problem setup to make the elements fully functional
237  // ----------------------------------------------------------------
238 
239  // Determine number of bulk elements in lower/upper fluids
240  const unsigned n_lower = Bulk_mesh_pt->nlower();
241  const unsigned n_upper = Bulk_mesh_pt->nupper();
242 
243  // Loop over bulk elements in lower fluid
244  for(unsigned e=0;e<n_lower;e++)
245  {
246  // Upcast from GeneralisedElement to the present element
247  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->
248  lower_layer_element_pt(e));
249 
250  // Set the Reynolds number
251  el_pt->re_pt() = &Global_Physical_Variables::Re;
252 
253  // Set the Womersley number
254  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
255 
256  // Set the product of the Reynolds number and the inverse of the
257  // Froude number
258  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
259 
260  // Set the direction of gravity
261  el_pt->g_pt() = &Global_Physical_Variables::G;
262 
263  } // End of loop over bulk elements in lower fluid
264 
265  // Loop over bulk elements in upper fluid
266  for(unsigned e=0;e<n_upper;e++)
267  {
268  // Upcast from GeneralisedElement to the present element
269  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->
270  upper_layer_element_pt(e));
271 
272  // Set the Reynolds number
273  el_pt->re_pt() = &Global_Physical_Variables::Re;
274 
275  // Set the Womersley number
276  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
277 
278  // Set the product of the Reynolds number and the inverse of the
279  // Froude number
280  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
281 
282  // Set the viscosity ratio
283  el_pt->viscosity_ratio_pt() = &Global_Physical_Variables::Viscosity_Ratio;
284 
285  // Set the density ratio
286  el_pt->density_ratio_pt() = &Global_Physical_Variables::Density_Ratio;
287 
288  // Set the direction of gravity
289  el_pt->g_pt() = &Global_Physical_Variables::G;
290 
291  } // End of loop over bulk elements in upper fluid
292 
293  // Set the pressure in the first element at 'node' 0 to 0.0
294  fix_pressure(0,0,0.0);
295 
296  // Determine number of 1D interface elements in mesh
297  unsigned n_interface_element = Surface_mesh_pt->nelement();
298 
299  // Loop over interface elements
300  for(unsigned e=0;e<n_interface_element;e++)
301  {
302  // Upcast from GeneralisedElement to the present element
306 
307  // Set the Strouhal number
309 
310  // Set the Capillary number
312 
313  } // End of loop over interface elements
314 
315  // Apply the boundary conditions
317 
318 
319  this->add_sub_mesh(Bulk_mesh_pt);
321 
322  this->build_global_mesh();
323 
324 
325  // Setup equation numbering scheme
326  cout << "Number of equations: " << assign_eqn_numbers() << std::endl;
327 
328 } // End of constructor
int i
Definition: BiCGSTAB_step_by_step.cpp:9
Definition: elements.h:1313
Definition: two_layer_spine_mesh.template.h:61
double St
Strouhal number.
Definition: elastic_two_layer_interface_axisym.cc:66

References oomph::Mesh::add_element_pt(), oomph::Problem::add_sub_mesh(), oomph::Problem::add_time_stepper_pt(), oomph::Problem::assign_eqn_numbers(), b, oomph::Problem::build_global_mesh(), InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt, Global_Physical_Variables::Ca, oomph::FluidInterfaceElement::ca_pt(), Global_Physical_Variables::Density_Ratio, e(), oomph::Mesh::element_pt(), InterfaceProblem< ELEMENT, TIMESTEPPER >::fix_pressure(), Global_Physical_Variables::G, i, n, oomph::Mesh::nelement(), Global_Physical_Variables::Re, Global_Physical_Variables::ReInvFr, Global_Physical_Variables::ReSt, InterfaceProblem< ELEMENT, TIMESTEPPER >::set_boundary_conditions(), Global_Physical_Variables::St, InterfaceProblem< ELEMENT, TIMESTEPPER >::Surface_mesh_pt, oomph::Problem::time_stepper_pt(), and Global_Physical_Variables::Viscosity_Ratio.

◆ ~InterfaceProblem() [4/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

108 {}

◆ InterfaceProblem() [5/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_r,
const unsigned n_y,
const unsigned n_theta,
const double r_min,
const double r_max,
const double l_y,
const double theta_max 
)

Problem constructor.

Constructor: Pass number of elements in x and y directions. Also lengths of the domain in x- and y-directions and the height of the layer

274  : R_max(r_max), L_y(l_y)
275 {
276  //this->linear_solver_pt() = new HSL_MA42;
277 
278  //static_cast<HSL_MA42*>(this->linear_solver_pt())->lenbuf_factor0() = 3.0;
279  //static_cast<HSL_MA42*>(this->linear_solver_pt())->lenbuf_factor1() = 3.0;
280  //static_cast<HSL_MA42*>(this->linear_solver_pt())->lenbuf_factor2() = 3.0;
281 
282  //Allocate the timestepper
283  add_time_stepper_pt(new TIMESTEPPER);
284 
285  //Now create the bulk mesh -- this should be your new annular mesh
287  (n_r,n_y,n_theta,r_min,r_max,l_y,0.0,theta_max,time_stepper_pt());
288 
289  //Set the documented node
290  Document_node_pt = Bulk_mesh_pt->boundary_node_pt(5,0);
291 
292 
293  //Create the surface mesh that will contain the interface elements
294  //First create storage, but with no elements or nodes
295  Surface_mesh_pt = new Mesh;
296 
297  // Loop over those elements adjacent to the free surface,
298  // which we shall choose to be the upper surface
299  for(unsigned e1=0;e1<n_y;e1++)
300  {
301  for(unsigned e2=0;e2<n_theta;e2++)
302  {
303  // Set a pointer to the bulk element we wish to our interface
304  // element to
305  FiniteElement* bulk_element_pt =
306  Bulk_mesh_pt->finite_element_pt(n_theta*n_y*(n_r-1) + e2 + e1*n_theta);
307 
308  // Create the interface element (on face 3 of the bulk element)
309  FiniteElement* interface_element_pt =
311 
312  // Add the interface element to the surface mesh
313  this->Surface_mesh_pt->add_element_pt(interface_element_pt);
314  }
315  }
316 
317  // Add the two sub-meshes to the problem
320 
321  // Combine all sub-meshes into a single mesh
323 
324  //Pin all nodes on the bottom
325  unsigned long n_boundary_node = Bulk_mesh_pt->nboundary_node(0);
326  for(unsigned long n=0;n<n_boundary_node;n++)
327  {
328  for(unsigned i=0;i<3;i++)
329  {
330  Bulk_mesh_pt->boundary_node_pt(0,n)->pin(i);
331  }
332  }
333 
334  //On the front and back (y=const) pin in y-direction
335  for(unsigned b=1;b<4;b+=2)
336  {
337  n_boundary_node = Bulk_mesh_pt->nboundary_node(b);
338  for(unsigned long n=0;n<n_boundary_node;n++)
339  {
340  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(1);
341  }
342  }
343 
344  //On sides pin in z-direction
345  for(unsigned b=4;b<5;b+=2)
346  {
347  n_boundary_node = Bulk_mesh_pt->nboundary_node(b);
348  for(unsigned long n=0;n<n_boundary_node;n++)
349  {
350  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(2);
351  }
352  }
353 
354  // pinning the top surface
355  for(unsigned b=2;b<3;b++)
356  {
357  n_boundary_node = Bulk_mesh_pt->nboundary_node(b);
358  for(unsigned long n=0;n<n_boundary_node;n++)
359  {
360  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(0);
361  }
362  }
363 
364  //Loop over the upper surface and set the surface concentration
365  {
366  unsigned b=5;
367  n_boundary_node = Bulk_mesh_pt->nboundary_node(b);
368  for(unsigned long n=0;n<n_boundary_node;n++)
369  {
370  Bulk_mesh_pt->boundary_node_pt(b,n)->set_value(3,1.0);;
371  }
372  }
373 
374 
375  //Create a Data object whose single value stores the
376  //external pressure
377  Data* external_pressure_data_pt = new Data(1);
378 
379  // Set and pin the external pressure to some random value
380  external_pressure_data_pt->set_value(0,1.31);
381  external_pressure_data_pt->pin(0);
382 
383  //Complete the problem setup to make the elements fully functional
384 
385  //Loop over the elements in the layer
386  unsigned n_bulk=Bulk_mesh_pt->nelement();
387  for(unsigned e=0;e<n_bulk;e++)
388  {
389  //Cast to a fluid element
390  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->element_pt(e));
391 
392  //Set the Reynolds number, etc
393  el_pt->re_pt() = &Global_Physical_Variables::Re;
394  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
395  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
396  el_pt->g_pt() = &Global_Physical_Variables::G;
397  }
398 
399  //Loop over 2D interface elements and set capillary number and
400  //the external pressure
401  unsigned long interface_element_pt_range =
403  for(unsigned e=0;e<interface_element_pt_range;e++)
404  {
405  //Cast to a interface element
409 
410  //Set the Capillary number
412 
413  //Pass the Data item that contains the single external pressure value
414  el_pt->set_external_pressure_data(external_pressure_data_pt);
415 
416  // Set the surface elasticity number
418 
419  // Set the surface peclect number
421 
422  // Set the surface peclect number multiplied by strouhal number
424 
425  }
426 
427  //Do equation numbering
428  cout << assign_eqn_numbers() << std::endl;
429 
430  //Now sort the elements to minimise frontwidth
431  std::sort(mesh_pt()->element_pt().begin(),
432  mesh_pt()->element_pt().end(),
433  ElementCmp());
434 
435 }
Function-type-object to perform comparison of elements.
Definition: elastic_bretherton.cc:317
Node * Document_node_pt
Pointer to a node for documentation purposes.
Definition: 3d_rayleigh_instability_surfactant.cc:261
double R_max
Axial lengths of domain.
Definition: 3d_rayleigh_instability_surfactant.cc:250
double L_y
Definition: 3d_rayleigh_instability_surfactant.cc:252
Definition: 3d_rayleigh_instability_surfactant.cc:121
Mesh *& mesh_pt()
Return a pointer to the global mesh.
Definition: problem.h:1280
Specialise to surface geometry.
Definition: surfactant_transport_elements.h:350
double *& peclet_s_pt()
Access function for pointer to the surface Peclet number.
Definition: surfactant_transport_elements.h:159
double *& beta_pt()
Access function for pointer to the Elasticity number.
Definition: surfactant_transport_elements.h:153
double *& peclet_strouhal_s_pt()
Access function for pointer to the surface Peclet x Strouhal number.
Definition: surfactant_transport_elements.h:165
static constexpr lastp1_t end
Definition: IndexedViewHelper.h:79
double r_min
Definition: two_d_biharmonic.cc:237
double r_max
Definition: two_d_biharmonic.cc:238
double Beta
Surface Elasticity number.
Definition: surfactant.cc:622
double Peclet_S
Surface Peclet number.
Definition: surfactant.cc:625
double Peclet_St_S
\shorT Sufrace Peclet number multiplied by Strouhal number
Definition: surfactant.cc:628

References oomph::Mesh::add_element_pt(), oomph::Problem::add_sub_mesh(), oomph::Problem::add_time_stepper_pt(), oomph::Problem::assign_eqn_numbers(), b, Global_Physical_Variables::Beta, oomph::SurfactantTransportInterfaceElement::beta_pt(), oomph::Problem::build_global_mesh(), InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt, Global_Physical_Variables::Ca, oomph::FluidInterfaceElement::ca_pt(), InterfaceProblem< ELEMENT, TIMESTEPPER >::Document_node_pt, e(), oomph::Mesh::element_pt(), Eigen::placeholders::end, Global_Physical_Variables::G, i, oomph::Problem::mesh_pt(), n, oomph::Mesh::nelement(), Global_Physical_Variables::Peclet_S, oomph::SurfactantTransportInterfaceElement::peclet_s_pt(), Global_Physical_Variables::Peclet_St_S, oomph::SurfactantTransportInterfaceElement::peclet_strouhal_s_pt(), oomph::Data::pin(), BiharmonicTestFunctions2::r_max, BiharmonicTestFunctions2::r_min, Global_Physical_Variables::Re, Global_Physical_Variables::ReInvFr, Global_Physical_Variables::ReSt, oomph::FluidInterfaceElement::set_external_pressure_data(), oomph::Data::set_value(), InterfaceProblem< ELEMENT, TIMESTEPPER >::Surface_mesh_pt, and oomph::Problem::time_stepper_pt().

◆ InterfaceProblem() [6/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_r,
const unsigned n_z,
const double l_z 
)

Constructor for single fluid interface problem.

Constructor: Pass the number of elements in radial and axial directions and the length of the domain in the z direction)

920 {
921  // Allocate the timestepper (this constructs the time object as well)
922  add_time_stepper_pt(new TIMESTEPPER(true));
923 
924  // Build and assign mesh (the "false" boolean flag tells the mesh
925  // constructor that the domain is not periodic in r)
927 
928  //Create "surface mesh" that will only contain the interface elements
929  Interface_mesh_pt = new Mesh;
930  {
931  // How many bulk elements are adjacent to boundary b?
932  // Boundary 1 is the outer boundary
933  // The boundaries are labelled
934  // 2
935  // 3 1
936  // 0
937 
938  unsigned n_element = Bulk_mesh_pt->nboundary_element(3);
939 
940  // Loop over the bulk elements adjacent to boundary b?
941  for(unsigned e=0;e<n_element;e++)
942  {
943  // Get pointer to the bulk element that is adjacent to boundary b
944  ELEMENT* bulk_elem_pt = dynamic_cast<ELEMENT*>(
945  Bulk_mesh_pt->boundary_element_pt(3,e));
946 
947  // Find the index of the face of element e along boundary b
948  int face_index = Bulk_mesh_pt->face_index_at_boundary(3,e);
949 
950  // Build the corresponding free surface element
953 
954  //Add the prescribed-flux element to the surface mesh
955  Interface_mesh_pt->add_element_pt(interface_element_pt);
956 
957  } //end of loop over bulk elements adjacent to boundary b
958  }
959 
960 
961  // Add the two sub meshes to the problem
964 
965  // Combine all submeshes into a single Mesh
967 
968 
969  // --------------------------------------------
970  // Set the boundary conditions for this problem
971  // --------------------------------------------
972 
973  //Pin all azimuthal velocities
974  {
975  unsigned n_node = this->Bulk_mesh_pt->nnode();
976  for(unsigned n=0;n<n_node;++n)
977  {
978  this->Bulk_mesh_pt->node_pt(n)->pin(2);
979  }
980  }
981 
982  // All nodes are free by default -- just pin the ones that have
983  // Dirichlet conditions here
984  unsigned ibound = 3;
985  unsigned n_node = Bulk_mesh_pt->nboundary_node(ibound);
986  for(unsigned n=0;n<n_node;n++)
987  {
988  Bulk_mesh_pt->boundary_node_pt(ibound,n)->set_value(3,1.0);
989  }
990 
991  // Determine number of mesh boundaries
992  const unsigned n_boundary = Bulk_mesh_pt->nboundary();
993 
994  // Loop over mesh boundaries
995  for(unsigned b=0;b<n_boundary;b++)
996  {
997  // Determine number of nodes on boundary b
998  const unsigned n_node = Bulk_mesh_pt->nboundary_node(b);
999 
1000  // Loop over nodes on boundary b
1001  for(unsigned n=0;n<n_node;n++)
1002  {
1003  // Pin azimuthal velocity on bounds
1004  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(2);
1005 
1006  // Pin velocity on the outer wall
1007  if(b==1)
1008  {
1009  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(0);
1010  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(1);
1011  }
1012  // Pin axial velocity on bottom and top walls (no penetration)
1013  if(b==0 || b==2)
1014  {
1015  Bulk_mesh_pt->boundary_node_pt(b,n)->pin(1);
1016  }
1017  } // End of loop over nodes on boundary b
1018  } // End of loop over mesh boundaries
1019 
1020  // ----------------------------------------------------------------
1021  // Complete the problem setup to make the elements fully functional
1022  // ----------------------------------------------------------------
1023 
1024  // Determine number of bulk elements in mesh
1025  const unsigned n_bulk = Bulk_mesh_pt->nelement();
1026 
1027  // Loop over the bulk elements
1028  for(unsigned e=0;e<n_bulk;e++)
1029  {
1030  // Upcast from GeneralisedElement to the present element
1031  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(Bulk_mesh_pt->element_pt(e));
1032 
1033  // Set the Reynolds number
1034  el_pt->re_pt() = &Global_Physical_Variables::Re;
1035 
1036  // Set the Womersley number
1037  el_pt->re_st_pt() = &Global_Physical_Variables::ReSt;
1038 
1039  // Set the product of the Reynolds number and the inverse of the
1040  // Froude number
1041  el_pt->re_invfr_pt() = &Global_Physical_Variables::ReInvFr;
1042 
1043  // Set the direction of gravity
1044  el_pt->g_pt() = &Global_Physical_Variables::G;
1045 
1046  } // End of loop over bulk elements
1047 
1048  // Create a Data object whose single value stores the external pressure
1049  Data* external_pressure_data_pt = new Data(1);
1050 
1051  // Pin and set the external pressure to some arbitrary value
1052  double p_ext = Global_Physical_Variables::P_ext;
1053 
1054  external_pressure_data_pt->pin(0);
1055  external_pressure_data_pt->set_value(0,p_ext);
1056 
1057  // Determine number of 1D interface elements in mesh
1058  const unsigned n_interface_element = Interface_mesh_pt->nelement();
1059 
1060  // Loop over the interface elements
1061  for(unsigned e=0;e<n_interface_element;e++)
1062  {
1063  // Upcast from GeneralisedElement to the present element
1067 
1068  // Set the Capillary number
1070 
1071  // Pass the Data item that contains the single external pressure value
1072  el_pt->set_external_pressure_data(external_pressure_data_pt);
1073 
1074  // Set the surface elasticity number
1076 
1077  // Set the surface peclect number
1079 
1080  // Set the surface peclect number multiplied by strouhal number
1082 
1083  } // End of loop over interface elements
1084 
1085  // Setup equation numbering scheme
1086  cout << "Number of equations: " << assign_eqn_numbers() << std::endl;
1087 
1088 } // End of constructor
Definition: rayleigh_instability_insoluble_surfactant.cc:711
Definition: rayleigh_instability_insoluble_surfactant.cc:117
double *& peclet_strouhal_s_pt()
Access function for pointer to the surface Peclet x Strouhal number.
Definition: rayleigh_instability_insoluble_surfactant.cc:573
double *& beta_pt()
Access function for pointer to the Elasticity number.
Definition: rayleigh_instability_insoluble_surfactant.cc:567
double *& peclet_s_pt()
Access function for pointer to the surface Peclet number.
Definition: rayleigh_instability_insoluble_surfactant.cc:570

References b, Global_Physical_Variables::Beta, oomph::SpineAxisymmetricMarangoniSurfactantFluidInterfaceElement< ELEMENT >::beta_pt(), oomph::Global_Physical_Variables::Ca, oomph::FluidInterfaceElement::ca_pt(), e(), Global_Physical_Variables::G, n, oomph::Mesh::nboundary_element(), Global_Physical_Variables::P_ext, Global_Physical_Variables::Peclet_S, oomph::SpineAxisymmetricMarangoniSurfactantFluidInterfaceElement< ELEMENT >::peclet_s_pt(), Global_Physical_Variables::Peclet_St_S, oomph::SpineAxisymmetricMarangoniSurfactantFluidInterfaceElement< ELEMENT >::peclet_strouhal_s_pt(), oomph::Data::pin(), oomph::Global_Physical_Variables::Re, oomph::Global_Physical_Variables::ReInvFr, Global_Physical_Variables::ReSt, oomph::FluidInterfaceElement::set_external_pressure_data(), and oomph::Data::set_value().

◆ ~InterfaceProblem() [5/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

761 {}

◆ InterfaceProblem() [7/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_r,
const unsigned n_z,
const double l_z 
)

Constructor: Pass the number of elements in radial and axial directions and the length of the domain in the z direction)

◆ ~InterfaceProblem() [6/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

424 {}

◆ InterfaceProblem() [8/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_x,
const unsigned n_y,
const double l_x,
const double h 
)

Constructor: Pass the number of elements and the lengths of the domain in the x and y directions (h is the height of the fluid layer i.e. the length of the domain in the y direction)

◆ ~InterfaceProblem() [7/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

104 {}

◆ InterfaceProblem() [9/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_x,
const unsigned n_y,
const double l_x,
const double h 
)

Constructor: Pass the number of elements and the lengths of the domain in the x and y directions (h is the height of the fluid layer i.e. the length of the domain in the y direction)

◆ ~InterfaceProblem() [8/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

95 {}

◆ InterfaceProblem() [10/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned Nx,
const unsigned Ny,
const unsigned Nz,
const double Lx,
const double Ly,
const double h 
)

Constructor: Pass number of elements in x and y directions. Also lengths of the domain in x- and y-directions and the height of the layer

◆ InterfaceProblem() [11/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( )

Constructor.

◆ ~InterfaceProblem() [9/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

182 {}

◆ InterfaceProblem() [12/12]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem ( const unsigned n_x,
const unsigned n_y1,
const unsigned n_y2,
const double l_x,
const double h1,
const double h2 
)

Constructor: Pass the number of elements and the width of the domain in the x direction. Also pass the number of elements in the y direction of the bottom (fluid 1) and top (fluid 2) layers, along with the heights of both layers.

◆ ~InterfaceProblem() [10/10]

template<class ELEMENT , class TIMESTEPPER >
InterfaceProblem< ELEMENT, TIMESTEPPER >::~InterfaceProblem ( )
inline

Destructor (empty)

105 {}

Member Function Documentation

◆ actions_after_adapt() [1/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_adapt
privatevirtual

Rebuild the mesh of interface elements after adapting the bulk mesh.

Reimplemented from oomph::Problem.

538 {
539  // Create the interface elements
541 
542  // Rebuild the Problem's global mesh from its various sub-meshes
544 
545  // Pin horizontal displacement of all nodes
546  const unsigned n_node = Bulk_mesh_pt->nnode();
547  for(unsigned n=0;n<n_node;n++) { Bulk_mesh_pt->node_pt(n)->pin_position(0); }
548 
549  // Unpin all fluid pressure dofs
550  RefineableAxisymmetricNavierStokesEquations::
551  unpin_all_pressure_dofs(Bulk_mesh_pt->element_pt());
552 
553  // Pin redudant fluid pressure dofs
554  RefineableAxisymmetricNavierStokesEquations::
555  pin_redundant_nodal_pressures(Bulk_mesh_pt->element_pt());
556 
557  // Now set the pressure in the first element at 'node' 0 to 0.0
558  fix_pressure(0,0,0.0);
559 
560  // Pin the redundant solid pressures (if any)
562  Bulk_mesh_pt->element_pt());
563 
564  // Reset the boundary conditions
566 
567 } // End of actions_after_adapt
void rebuild_global_mesh()
Definition: problem.cc:1533

References n.

◆ actions_after_adapt() [2/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_adapt ( )
privatevirtual

Rebuild the mesh of interface elements after adapting the bulk mesh.

Reimplemented from oomph::Problem.

◆ actions_after_newton_solve() [1/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlinevirtual

Update after solve can remain empty, because the update is performed automatically after every Newton step.

Reimplemented from oomph::Problem.

116 {}

◆ actions_after_newton_solve() [2/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlinevirtual

Update after solve can remain empty, because the update is performed automatically after every Newton step.

Reimplemented from oomph::Problem.

109 {}

◆ actions_after_newton_solve() [3/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlineprivatevirtual

No actions required after solve step.

Reimplemented from oomph::Problem.

202 {}

◆ actions_after_newton_solve() [4/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlineprivatevirtual

Update after solve can remain empty, because the update is performed automatically after every Newton step.

Reimplemented from oomph::Problem.

138 {}

◆ actions_after_newton_solve() [5/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlineprivatevirtual

No actions required after solve step.

Reimplemented from oomph::Problem.

124 {}

◆ actions_after_newton_solve() [6/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlineprivatevirtual

Update after solve can remain empty, because the update is performed automatically after every Newton step.

Reimplemented from oomph::Problem.

131 {}

◆ actions_after_newton_solve() [7/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlineprivatevirtual

No actions required after solve step.

Reimplemented from oomph::Problem.

202 {}

◆ actions_after_newton_solve() [8/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_after_newton_solve ( )
inlineprivatevirtual

Update after solve can remain empty, because the update is performed automatically after every Newton step.

Reimplemented from oomph::Problem.

135 {}

◆ actions_before_adapt() [1/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_adapt
privatevirtual

Strip off the interface elements before adapting the bulk mesh.

Reimplemented from oomph::Problem.

522 {
523  // Delete the interface elements and wipe the surface mesh
525 
526  // Rebuild the Problem's global mesh from its various sub-meshes
528 
529 } // End of actions_before_adapt
void delete_interface_elements()
Delete the 1d interface elements.
Definition: elastic_two_layer_interface_axisym.cc:639

◆ actions_before_adapt() [2/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_adapt ( )
privatevirtual

Strip off the interface elements before adapting the bulk mesh.

Reimplemented from oomph::Problem.

◆ actions_before_implicit_timestep() [1/3]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_implicit_timestep ( )
inlineprivatevirtual

Actions before the timestep: For maximum stability, reset the current nodal positions to be the "stress-free" ones.

Reimplemented from oomph::Problem.

207  {
208  Bulk_mesh_pt->set_lagrangian_nodal_coordinates();
209  }

◆ actions_before_implicit_timestep() [2/3]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_implicit_timestep ( )
inlineprivatevirtual

Actions before the timestep: For maximum stability, reset the current nodal positions to be the "stress-free" ones.

Reimplemented from oomph::Problem.

129  {
130  Bulk_mesh_pt->set_lagrangian_nodal_coordinates();
131  }

◆ actions_before_implicit_timestep() [3/3]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_implicit_timestep ( )
inlineprivatevirtual

Actions before the timestep: For maximum stability, reset the current nodal positions to be the "stress-free" ones.

Reimplemented from oomph::Problem.

207  {
208  Bulk_mesh_pt->set_lagrangian_nodal_coordinates();
209  }

◆ actions_before_newton_convergence_check() [1/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlinevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them here.

Reimplemented from oomph::Problem.

107  {
108  Bulk_mesh_pt->node_update();
109  }

◆ actions_before_newton_convergence_check() [2/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlinevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them here.

Reimplemented from oomph::Problem.

100  {
101  Bulk_mesh_pt->node_update();
102  }

◆ actions_before_newton_convergence_check() [3/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlineprivatevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them here.

Reimplemented from oomph::Problem.

129  {
130  Bulk_mesh_pt->node_update();
131  }

◆ actions_before_newton_convergence_check() [4/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlinevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them

Reimplemented from oomph::Problem.

213 {Bulk_mesh_pt->node_update();}

◆ actions_before_newton_convergence_check() [5/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlinevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them here.

Reimplemented from oomph::Problem.

768  {
769  Bulk_mesh_pt->node_update();
770  }

◆ actions_before_newton_convergence_check() [6/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlineprivatevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them here.

Reimplemented from oomph::Problem.

122  {
123  Bulk_mesh_pt->node_update();
124  }

References oomph::SpineMesh::node_update().

◆ actions_before_newton_convergence_check() [7/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlinevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them

Reimplemented from oomph::Problem.

102 {Bulk_mesh_pt->node_update();}

◆ actions_before_newton_convergence_check() [8/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_convergence_check ( )
inlineprivatevirtual

Spine heights/lengths are unknowns in the problem so their values get corrected during each Newton step. However, changing their value does not automatically change the nodal positions, so we need to update all of them here.

Reimplemented from oomph::Problem.

126  {
127  Bulk_mesh_pt->node_update();
128  }

◆ actions_before_newton_solve() [1/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlinevirtual

Any actions that are to be performed before a complete Newton solve (e.g. adjust boundary conditions). CAREFUL: This step should (and if the FD-based LinearSolver FD_LU is used, must) only update values that are pinned!

Reimplemented from oomph::Problem.

112 {}

◆ actions_before_newton_solve() [2/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlinevirtual

Any actions that are to be performed before a complete Newton solve (e.g. adjust boundary conditions). CAREFUL: This step should (and if the FD-based LinearSolver FD_LU is used, must) only update values that are pinned!

Reimplemented from oomph::Problem.

105 {}

◆ actions_before_newton_solve() [3/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlineprivatevirtual

No actions required before solve step.

Reimplemented from oomph::Problem.

199 {}

◆ actions_before_newton_solve() [4/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlineprivatevirtual

No actions required before solve step.

Reimplemented from oomph::Problem.

134 {}

◆ actions_before_newton_solve() [5/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlineprivatevirtual

No actions required before solve step.

Reimplemented from oomph::Problem.

121 {}

◆ actions_before_newton_solve() [6/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlineprivatevirtual

No actions required before solve step.

Reimplemented from oomph::Problem.

127 {}

◆ actions_before_newton_solve() [7/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlineprivatevirtual

No actions required before solve step.

Reimplemented from oomph::Problem.

199 {}

◆ actions_before_newton_solve() [8/8]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::actions_before_newton_solve ( )
inlineprivatevirtual

No actions required before solve step.

Reimplemented from oomph::Problem.

131 {}

◆ compute_total_mass() [1/3]

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::compute_total_mass ( )
inline

Compute the total mass.

223  {
224  double mass = 0.0;
225 
226  // Determine number of 1D interface elements in mesh
227  const unsigned n_interface_element = Surface_mesh_pt->nelement();
228 
229  // Loop over the interface elements
230  for(unsigned e=0;e<n_interface_element;e++)
231  {
232  // Upcast from GeneralisedElement to the present element
236 
237  mass += el_pt->integrate_c();
238  }
239  return mass;
240  }
double integrate_c()
Compute the concentration intergated over the surface area.
Definition: surfactant_transport_elements.cc:358

References e(), and oomph::SurfactantTransportInterfaceElement::integrate_c().

◆ compute_total_mass() [2/3]

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::compute_total_mass ( )
inline

Compute the total mass of the insoluble surfactant.

856  {
857  //Initialise to zero
858  double mass = 0.0;
859 
860  // Determine number of 1D interface elements in mesh
861  const unsigned n_interface_element = Interface_mesh_pt->nelement();
862 
863  // Loop over the interface elements
864  for(unsigned e=0;e<n_interface_element;e++)
865  {
866  // Upcast from GeneralisedElement to the present element
870  //Add contribution from each element
871  mass += el_pt->integrate_c();
872  }
873  return mass;
874  } // End of compute_total_mass
double integrate_c() const
Compute the concentration intergated over the area.
Definition: rayleigh_instability_insoluble_surfactant.cc:627

References e(), oomph::Mesh::element_pt(), oomph::SpineAxisymmetricMarangoniSurfactantFluidInterfaceElement< ELEMENT >::integrate_c(), and oomph::Mesh::nelement().

◆ compute_total_mass() [3/3]

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::compute_total_mass ( )
inline

Compute the total mass of the insoluble surfactant.

518  {
519  //Initialise to zero
520  double mass = 0.0;
521 
522  // Determine number of 1D interface elements in mesh
523  const unsigned n_interface_element = Interface_mesh_pt->nelement();
524 
525  // Loop over the interface elements
526  for(unsigned e=0;e<n_interface_element;e++)
527  {
528  // Upcast from GeneralisedElement to the present element
532  //Add contribution from each element
533  mass += el_pt->integrate_c();
534  }
535  return mass;
536  } // End of compute_total_mass
Specialise to the Axisymmetric geometry.
Definition: surfactant_transport_elements.h:312

References e(), oomph::Mesh::element_pt(), oomph::SurfactantTransportInterfaceElement::integrate_c(), and oomph::Mesh::nelement().

◆ create_interface_elements() [1/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::create_interface_elements
private

Create the 1d interface elements.

Create interface elements between the two fluids in the mesh pointed to by Bulk_mesh_pt and add the elements to the Mesh object pointed to by Surface_mesh_pt.

578 {
579  // Determine number of bulk elements adjacent to interface (boundary 4)
580  const unsigned n_element = this->Bulk_mesh_pt->nboundary_element(4);
581 
582  // Loop over those elements adjacent to the interface
583  for(unsigned e=0;e<n_element;e++)
584  {
585  // Get pointer to the bulk element that is adjacent to the interface
586  ELEMENT* bulk_elem_pt = dynamic_cast<ELEMENT*>(
587  this->Bulk_mesh_pt->boundary_element_pt(4,e));
588 
589  // We only want to attach interface elements to the bulk elements
590  // which are BELOW the interface, and so we filter out those above by
591  // referring to the viscosity_ratio_pt
592  if(bulk_elem_pt->viscosity_ratio_pt()
594  {
595  // Find index of the face of element e that corresponds to the interface
596  const int face_index = this->Bulk_mesh_pt->face_index_at_boundary(4,e);
597 
598  // Create the interface element
599  FiniteElement* interface_element_element_pt =
601  face_index);
602 
603  // Add the interface element to the surface mesh
604  this->Surface_mesh_pt->add_element_pt(interface_element_element_pt);
605  }
606  }
607 
608  // --------------------------------------------------------
609  // Complete the setup to make the elements fully functional
610  // --------------------------------------------------------
611 
612  // Determine number of 1D interface elements in mesh
613  const unsigned n_interface_element = this->Surface_mesh_pt->nelement();
614 
615  // Loop over the interface elements
616  for(unsigned e=0;e<n_interface_element;e++)
617  {
618  // Upcast from GeneralisedElement to the present element
622 
623  // Set the Strouhal number
625 
626  // Set the Capillary number
628 
629  } // End of loop over interface elements
630 
631 } // End of create_interface_elements()
Specialise the Elastic update case to axisymmetric equations.
Definition: specific_node_update_interface_elements.h:1257
double *& st_pt()
The pointer to the Strouhal number.
Definition: interface_elements.h:504

References Global_Physical_Variables::Ca, oomph::FluidInterfaceElement::ca_pt(), e(), Global_Physical_Variables::St, oomph::FluidInterfaceElement::st_pt(), and Global_Physical_Variables::Viscosity_Ratio.

◆ create_interface_elements() [2/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::create_interface_elements ( )
private

Create the 1d interface elements.

◆ deform_free_surface() [1/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon)
inlineprivate

Deform the mesh/free surface to a prescribed function.

161  {
162  // Determine number of spines in mesh
163  const unsigned n_spine = Bulk_mesh_pt->nspine();
164 
165  // Loop over spines in mesh
166  for(unsigned i=0;i<n_spine;i++)
167  {
168 
169  // Determine z coordinate of spine
170  double z_value = Bulk_mesh_pt->boundary_node_pt(1,i)->x(1);
171 
172  // Set spine height
173  Bulk_mesh_pt->spine_pt(i)->height() = Lr + epsilon*cos((z_value/Height)*2.0*MathematicalConstants::Pi);
174 
175  } // End of loop over spines
176 
177  // Update nodes in bulk mesh
178  Bulk_mesh_pt->node_update();
179 
180  } // End of deform_free_surface
AnnoyingScalar cos(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:136
double Pi
Definition: two_d_biharmonic.cc:235
double epsilon
Definition: osc_ring_sarah_asymptotics.h:43

References oomph::Mesh::boundary_node_pt(), cos(), oomph::SarahBL::epsilon, GlobalParameters::Height, oomph::Spine::height(), i, Global_Parameters::Lr, oomph::SpineMesh::node_update(), oomph::SpineMesh::nspine(), BiharmonicTestFunctions2::Pi, oomph::SpineMesh::spine_pt(), and oomph::Node::x().

◆ deform_free_surface() [2/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon)
inlineprivate

Deform the mesh/free surface to a prescribed function.

881  {
882  // Determine number of spines in mesh
883  const unsigned n_spine = Bulk_mesh_pt->nspine();
884 
885  // Loop over spines in mesh
886  for(unsigned i=0;i<n_spine;i++)
887  {
888 
889  // Determine z coordinate of spine
890  double z_value = Bulk_mesh_pt->boundary_node_pt(1,i)->x(1);
891 
892  // Set spine height
893  Bulk_mesh_pt->spine_pt(i)->height() =
896 
897  } // End of loop over spines
898 
899  // Update nodes in bulk mesh
900  Bulk_mesh_pt->node_update();
901 
902  } // End of deform_free_surface
double Film_Thickness
Definition: 3d_rayleigh_instability_surfactant.cc:54
double Alpha
Scaling factor for wall thickness (to be used in an exercise)
Definition: unsteady_ring.cc:73

References Global_Physical_Variables::Alpha, oomph::Mesh::boundary_node_pt(), cos(), oomph::SarahBL::epsilon, Global_Physical_Variables::Film_Thickness, oomph::Spine::height(), i, oomph::SpineMesh::node_update(), oomph::SpineMesh::nspine(), oomph::SpineMesh::spine_pt(), and oomph::Node::x().

◆ deform_free_surface() [3/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon)
inlineprivate

Deform the mesh/free surface to a prescribed function.

543  {
544  //Loop over all nodes in the mesh
545  const unsigned n_node = Bulk_mesh_pt->nnode();
546  for(unsigned n=0;n<n_node;n++)
547  {
548  //Find out the r value
549  double r = Bulk_mesh_pt->node_pt(n)->x(0);
550  // Determine z coordinate of spine
551  double z_value = Bulk_mesh_pt->node_pt(n)->x(1);
552 
553  //Set the thickess of the filme
554  double thickness =
557 
558  //Now scale the position accordingly
559  Bulk_mesh_pt->node_pt(n)->x(0) = 1.0 - (1.0-r)*thickness;
560  }
561 
562  //Reset the lagrangian coordinates
563  Bulk_mesh_pt->set_lagrangian_nodal_coordinates();
564  } // End of deform_free_surface
r
Definition: UniformPSDSelfTest.py:20

References Global_Physical_Variables::Alpha, cos(), oomph::SarahBL::epsilon, Global_Physical_Variables::Film_Thickness, n, oomph::Mesh::nnode(), oomph::SolidMesh::node_pt(), UniformPSDSelfTest::r, oomph::SolidMesh::set_lagrangian_nodal_coordinates(), and oomph::Node::x().

◆ deform_free_surface() [4/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon,
const double k 
)
inlineprivate

Deform the mesh/free surface to a prescribed function.

153  {
154  // Initialise Bessel functions (only need the first!)
155  double j0, j1, y0, y1, j0p, j1p, y0p, y1p;
156 
157  // Determine number of spines in mesh
158  const unsigned n_spine = Bulk_mesh_pt->nspine();
159 
160  // Loop over spines in mesh
161  for(unsigned i=0;i<n_spine;i++)
162  {
163 
164  // Determine r coordinate of spine
165  double r_value = Bulk_mesh_pt->boundary_node_pt(0,i)->x(0);
166 
167  // Get Bessel functions J_0(kr), J_1(kr), Y_0(kr), Y_1(kr)
168  // and their derivatives
169  CRBond_Bessel::bessjy01a(k*r_value,j0,j1,y0,y1,j0p,j1p,y0p,y1p);
170 
171  // Set spine height
172  Bulk_mesh_pt->spine_pt(i)->height() = 1.0 + epsilon*j0;
173 
174  } // End of loop over spines
175 
176  // Update nodes in bulk mesh
177  Bulk_mesh_pt->node_update();
178 
179  } // End of deform_free_surface
char char char int int * k
Definition: level2_impl.h:374
int bessjy01a(double x, double &j0, double &j1, double &y0, double &y1, double &j0p, double &j1p, double &y0p, double &y1p)
Definition: crbond_bessel.cc:551

References CRBond_Bessel::bessjy01a(), oomph::Mesh::boundary_node_pt(), oomph::SarahBL::epsilon, oomph::Spine::height(), i, k, oomph::SpineMesh::node_update(), oomph::SpineMesh::nspine(), oomph::SpineMesh::spine_pt(), and oomph::Node::x().

◆ deform_free_surface() [5/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon,
const double k 
)
private

Deform the mesh/free surface to a prescribed function.

◆ deform_free_surface() [6/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon,
const double k 
)
private

Deform the mesh/free surface to a prescribed function.

◆ deform_free_surface() [7/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon,
const unsigned n_periods 
)
private

Deform the mesh/free surface to a prescribed function.

405 {
406  // Determine number of nodes in the "bulk" mesh
407  const unsigned n_node = Bulk_mesh_pt->nnode();
408 
409  // Loop over all nodes in mesh
410  for(unsigned n=0;n<n_node;n++)
411  {
412  // Determine eulerian position of node
413  const double current_x_pos = Bulk_mesh_pt->node_pt(n)->x(0);
414  const double current_y_pos = Bulk_mesh_pt->node_pt(n)->x(1);
415 
416  // Determine new vertical position of node
417  const double new_y_pos = current_y_pos
418  + (1.0-fabs(1.0-current_y_pos))*epsilon
419  *(cos(2.0*n_periods*MathematicalConstants::Pi*current_x_pos/Lx));
420 
421  // Set new position
422  Bulk_mesh_pt->node_pt(n)->x(1) = new_y_pos;
423  }
424 } // End of deform_free_surface
double Lx
Width of domain.
Definition: elastic_single_layer.cc:146
Real fabs(const Real &a)
Definition: boostmultiprec.cpp:117

References cos(), oomph::SarahBL::epsilon, boost::multiprecision::fabs(), Global_Parameters::Lx, n, and BiharmonicTestFunctions2::Pi.

◆ deform_free_surface() [8/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon,
const unsigned n_periods 
)
private

Deform the mesh/free surface to a prescribed function.

◆ deform_free_surface() [9/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon,
const unsigned n_periods 
)
private

Deform the mesh/free surface to a prescribed function.

◆ deform_free_surface() [10/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::deform_free_surface ( const double epsilon,
const unsigned n_periods 
)
private

Deform the mesh/free surface to a prescribed function.

◆ delete_interface_elements() [1/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::delete_interface_elements
private

Delete the 1d interface elements.

Delete the interface elements and wipe the surface mesh.

640 {
641  // Determine number of interface elements
642  const unsigned n_interface_element = Surface_mesh_pt->nelement();
643 
644  // Loop over interface elements and delete
645  for(unsigned e=0;e<n_interface_element;e++)
646  {
647  delete Surface_mesh_pt->element_pt(e);
648  }
649 
650  // Wipe the mesh
652 
653 } // End of delete_interface_elements
void flush_element_and_node_storage()
Definition: mesh.h:407

References e().

◆ delete_interface_elements() [2/2]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::delete_interface_elements ( )
private

Delete the 1d interface elements.

◆ doc_solution() [1/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

Document the solution.

355 {
356 
357  // Output the time
358  double t= time_pt()->time();
359  cout << "Time is now " << t << std::endl;
360 
361  // Document in trace file
362  Trace_file << time_pt()->time() << " "
363  << Bulk_mesh_pt->spine_pt(0)->height() << std::endl;
364 
365  ofstream some_file;
366  char filename[100];
367 
368  // Set number of plot points (in each coordinate direction)
369  const unsigned npts = 5;
370 
371  // Open solution output file
372  sprintf(filename,"%s/soln%i.dat",doc_info.directory().c_str(),
373  doc_info.number());
374  some_file.open(filename);
375 
376  // Output solution to file
377  Bulk_mesh_pt->output(some_file,npts);
378  Interface_mesh_pt->output(some_file,npts);
379 
380  // Write file as a tecplot text object...
381  some_file << "TEXT X=2.5,Y=93.6,F=HELV,HU=POINT,C=BLUE,H=26,T=\"time = "
382  << time_pt()->time() << "\"";
383  // ...and draw a horizontal line whose length is proportional
384  // to the elapsed time
385  some_file << "GEOMETRY X=2.5,Y=98,T=LINE,C=BLUE,LT=0.4" << std::endl;
386  some_file << "1" << std::endl;
387  some_file << "2" << std::endl;
388  some_file << " 0 0" << std::endl;
389  some_file << time_pt()->time()*20.0 << " 0" << std::endl;
390 
391  // Close solution output file
392  some_file.close();
393 
394  // Output solution to file in paraview format
395  string file_name="soln"+StringConversion::to_string(doc_info.number())
396  +".vtu";
397  sprintf(filename,"%s/%s",doc_info.directory().c_str(),file_name.c_str());
398  some_file.open(filename);
399  Bulk_mesh_pt->output_paraview(some_file,npts);
400  some_file.close();
401 
402  // Write pvd information
404  file_name,t);
405 
406 } // End of doc_solution
ofstream Trace_file
Trace file.
Definition: rayleigh_instability.cc:183
std::string directory() const
Output directory.
Definition: oomph_utilities.h:524
unsigned & number()
Number used (e.g.) for labeling output files.
Definition: oomph_utilities.h:554
void output(std::ostream &outfile)
Output for all elements.
Definition: mesh.cc:2027
Time *& time_pt()
Return a pointer to the global time object.
Definition: problem.h:1504
double & time()
Return the current value of the continuous time.
Definition: timesteppers.h:123
ofstream Pvd_file
Definition: rayleigh_instability.cc:76
string filename
Definition: MergeRestartFiles.py:39
string file_name
Definition: Particles2023AnalysisHung.py:321
void write_pvd_information(std::ofstream &pvd_file, const std::string &output_filename, const double &time)
Definition: mesh.cc:9624
std::string to_string(T object, unsigned float_precision=8)
Definition: oomph_utilities.h:189
t
Definition: plotPSD.py:36

References oomph::DocInfo::directory(), Particles2023AnalysisHung::file_name, MergeRestartFiles::filename, oomph::DocInfo::number(), Global_Physical_Variables::Pvd_file, plotPSD::t, oomph::StringConversion::to_string(), oomph::Problem_Parameter::Trace_file, and oomph::ParaviewHelper::write_pvd_information().

◆ doc_solution() [2/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ doc_solution() [3/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Document the solution.

◆ doc_solution() [4/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ doc_solution() [5/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ doc_solution() [6/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ doc_solution() [7/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ doc_solution() [8/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Document the solution.

◆ doc_solution() [9/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ doc_solution() [10/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ doc_solution() [11/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Document the solution.

◆ doc_solution() [12/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::doc_solution ( DocInfo doc_info)

Doc the solution.

◆ fix_pressure() [1/4]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::fix_pressure ( const unsigned e,
const unsigned pdof,
const double pvalue 
)
inlineprivate

Fix pressure in element e at pressure dof pdof and set to pvalue.

230  {
231  // Fix the pressure at that element
232  dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(e))->
233  fix_pressure(pdof,pvalue);
234  }

References e().

Referenced by InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem().

◆ fix_pressure() [2/4]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::fix_pressure ( const unsigned e,
const unsigned pdof,
const double pvalue 
)
inlineprivate

Fix pressure in element e at pressure dof pdof and set to pvalue.

147  {
148  // Fix the pressure at that element
149  dynamic_cast<ELEMENT*>(Bulk_mesh_pt->element_pt(e))->
150  fix_pressure(pdof,pvalue);
151  }

References e().

◆ fix_pressure() [3/4]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::fix_pressure ( const unsigned e,
const unsigned pdof,
const double pvalue 
)
inlineprivate

Fix pressure in element e at pressure dof pdof and set to pvalue.

230  {
231  // Fix the pressure at that element
232  dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(e))->
233  fix_pressure(pdof,pvalue);
234  }

References e().

◆ fix_pressure() [4/4]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::fix_pressure ( const unsigned e,
const unsigned pdof,
const double pvalue 
)
inlineprivate

Fix pressure in element e at pressure dof pdof and set to pvalue.

144  {
145  // Fix the pressure at that element
146  dynamic_cast<ELEMENT*>(Bulk_mesh_pt->element_pt(e))->
147  fix_pressure(pdof,pvalue);
148  }

References e().

◆ global_temporal_error_norm() [1/2]

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::global_temporal_error_norm ( )
inlinevirtual

The global temporal error norm, based on the movement of the nodes in the radial direction only (because that's the only direction in which they move!)

Reimplemented from oomph::Problem.

802  {
803  //Temp
804  double global_error = 0.0;
805 
806  //Find out how many nodes there are in the problem
807  const unsigned n_node = Bulk_mesh_pt->nnode();
808 
809  //Loop over the nodes and calculate the errors in the positions
810  for(unsigned n=0;n<n_node;n++)
811  {
812  //Set the dimensions to be restricted to the radial direction only
813  const unsigned n_dim = 1;
814  //Set the position error to zero
815  double node_position_error = 0.0;
816  //Loop over the dimensions
817  for(unsigned i=0;i<n_dim;i++)
818  {
819  //Get position error
820  double error =
821  Bulk_mesh_pt->node_pt(n)->position_time_stepper_pt()->
822  temporal_error_in_position(Bulk_mesh_pt->node_pt(n),i);
823 
824  //Add the square of the individual error to the position error
825  node_position_error += error*error;
826  }
827 
828  //Divide the position error by the number of dimensions
829  node_position_error /= n_dim;
830  //Now add to the global error
831  global_error += node_position_error;
832  }
833 
834  //Now the global error must be divided by the number of nodes
835  global_error /= n_node;
836 
837  //Return the square root of the errr
838  return sqrt(global_error);
839  }
AnnoyingScalar sqrt(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:134
int error
Definition: calibrate.py:297

References calibrate::error, i, n, and sqrt().

◆ global_temporal_error_norm() [2/2]

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::global_temporal_error_norm ( )
inlinevirtual

The global temporal error norm, based on the movement of the nodes in the radial direction only (because that's the only direction in which they move!)

Reimplemented from oomph::Problem.

458  {
459  //Temp
460  double global_error = 0.0;
461 
462  //Find out how many nodes there are in the problem
463  const unsigned n_node = Bulk_mesh_pt->nnode();
464 
465  //Loop over the nodes and calculate the errors in the positions
466  for(unsigned n=0;n<n_node;n++)
467  {
468  //Set the dimensions to be restricted to the radial direction only
469  const unsigned n_dim = 1;
470  //Set the position error to zero
471  double node_position_error = 0.0;
472  //Loop over the dimensions
473  for(unsigned i=0;i<n_dim;i++)
474  {
475  //Get position error
476  double error =
477  Bulk_mesh_pt->node_pt(n)->position_time_stepper_pt()->
478  temporal_error_in_position(Bulk_mesh_pt->node_pt(n),i);
479 
480  //Add the square of the individual error to the position error
481  node_position_error += error*error;
482  }
483 
484  //Divide the position error by the number of dimensions
485  node_position_error /= n_dim;
486  //Now add to the global error
487  global_error += node_position_error;
488  }
489 
490  //Now the global error must be divided by the number of nodes
491  global_error /= n_node;
492 
493  //Return the square root of the errr
494  return sqrt(global_error);
495  }

References calibrate::error, i, n, and sqrt().

◆ set_boundary_conditions() [1/6]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_boundary_conditions

Set boundary conditions.

Set boundary conditions: Set all velocity components to zero on the top and bottom (solid) walls and the radial and azimuthal components only to zero on the side boundaries

Set boundary conditions: Set both velocity components to zero on the bottom (solid) wall and the horizontal component only to zero on the side (periodic) boundaries

Set boundary conditions: Set both velocity components to zero on the top and bottom (solid) walls and the horizontal component only to zero on the side (periodic) boundaries

485 {
486  // Determine number of mesh boundaries
487  const unsigned n_boundary = Bulk_mesh_pt->nboundary();
488 
489  // Loop over mesh boundaries
490  for(unsigned b=0;b<n_boundary;b++)
491  {
492  // Determine number of nodes on boundary b
493  const unsigned n_node = Bulk_mesh_pt->nboundary_node(b);
494 
495  // Loop over nodes on boundary b
496  for(unsigned n=0;n<n_node;n++)
497  {
498  // Set radial component of the velocity to zero on all boundaries
499  // other than the interface (b=4)
500  if(b!=4) { Bulk_mesh_pt->boundary_node_pt(b,n)->set_value(0,0.0); }
501 
502  // Set azimuthal component of the velocity to zero on all boundaries
503  // other than the interface (b=4)
504  if(b!=4) { Bulk_mesh_pt->boundary_node_pt(b,n)->set_value(2,0.0); }
505 
506  // Set axial component of the velocity to zero on solid boundaries
507  if(b==0 || b==2)
508  {
509  Bulk_mesh_pt->boundary_node_pt(b,n)->set_value(1,0.0);
510  }
511  }
512  }
513 } // End of set_boundary_conditions

References b, and n.

Referenced by InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem().

◆ set_boundary_conditions() [2/6]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_boundary_conditions ( )

Set boundary conditions.

◆ set_boundary_conditions() [3/6]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_boundary_conditions ( )

Set boundary conditions.

◆ set_boundary_conditions() [4/6]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_boundary_conditions ( )

Set boundary conditions.

◆ set_boundary_conditions() [5/6]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_boundary_conditions ( )

Set boundary conditions.

◆ set_boundary_conditions() [6/6]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_boundary_conditions ( )

Set boundary conditions.

◆ set_initial_condition() [1/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition
inlinevirtual

Set initial conditions: Set all nodal velocities to zero and initialise the previous velocities to correspond to an impulsive start

Set initial conditions: Set all nodal velocities to zero and initialise the previous velocities and nodal positions to correspond to an impulsive start

Reimplemented from oomph::Problem.

123  {
124  // Determine number of nodes in mesh
125  const unsigned n_node = Bulk_mesh_pt->nnode();
126 
127  // Loop over all nodes in mesh
128  for(unsigned n=0;n<n_node;n++)
129  {
130  // Loop over the three velocity components
131  for(unsigned i=0;i<3;i++)
132  {
133  // Set velocity component i of node n to zero
134  Bulk_mesh_pt->node_pt(n)->set_value(i,0.0);
135  }
136  }
137 
138  // Initialise the previous velocity values for timestepping
139  // corresponding to an impulsive start
141 
142  } // End of set_initial_condition
void assign_initial_values_impulsive()
Definition: problem.cc:11499

References i, and n.

◆ set_initial_condition() [2/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
inlinevirtual

Set initial conditions: Set all nodal velocities to zero and initialise the previous velocities to correspond to an impulsive start

Reimplemented from oomph::Problem.

115  {
116  // Determine number of nodes in mesh
117  const unsigned n_node = Bulk_mesh_pt->nnode();
118 
119  // Loop over all nodes in mesh
120  for(unsigned n=0;n<n_node;n++)
121  {
122  // Loop over the three velocity components
123  for(unsigned i=0;i<3;i++)
124  {
125  // Set velocity component i of node n to zero
126  Bulk_mesh_pt->node_pt(n)->set_value(i,0.0);
127  }
128  }
129 
130  // Initialise the previous velocity values for timestepping
131  // corresponding to an impulsive start
133 
134  } // End of set_initial_condition

References i, and n.

◆ set_initial_condition() [3/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
virtual

Set initial conditions.

Reimplemented from oomph::Problem.

◆ set_initial_condition() [4/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
virtual

Set initial conditions.

Reimplemented from oomph::Problem.

◆ set_initial_condition() [5/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
inlinevirtual

Set initial conditions: Set all nodal velocities to zero and initialise the previous velocities to correspond to an impulsive start

Reimplemented from oomph::Problem.

776  {
777  // Determine number of nodes in mesh
778  const unsigned n_node = Bulk_mesh_pt->nnode();
779 
780  // Loop over all nodes in mesh
781  for(unsigned n=0;n<n_node;n++)
782  {
783  // Loop over the three velocity components
784  for(unsigned i=0;i<3;i++)
785  {
786  // Set velocity component i of node n to zero
787  Bulk_mesh_pt->node_pt(n)->set_value(i,0.0);
788  }
789  }
790 
791  // Initialise the previous velocity values for timestepping
792  // corresponding to an impulsive start
794 
795  } // End of set_initial_condition

References i, and n.

◆ set_initial_condition() [6/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
inlinevirtual

Set initial conditions: Set all nodal velocities to zero and initialise the previous velocities to correspond to an impulsive start

Reimplemented from oomph::Problem.

430  {
431  // Determine number of nodes in mesh
432  const unsigned n_node = Bulk_mesh_pt->nnode();
433 
434  // Loop over all nodes in mesh
435  for(unsigned n=0;n<n_node;n++)
436  {
437  // Loop over the three velocity components
438  for(unsigned i=0;i<3;i++)
439  {
440  // Set velocity component i of node n to zero
441  Bulk_mesh_pt->node_pt(n)->set_value(i,0.0);
442  }
443  //Set the bulk concentration to be 1
444  Bulk_mesh_pt->node_pt(n)->set_value(3,1.0);
445  }
446 
447  // Initialise the previous velocity values for timestepping
448  // corresponding to an impulsive start
450 
451  } // End of set_initial_condition

References i, and n.

◆ set_initial_condition() [7/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
virtual

Set initial conditions.

Reimplemented from oomph::Problem.

◆ set_initial_condition() [8/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
virtual

Set initial conditions.

Reimplemented from oomph::Problem.

◆ set_initial_condition() [9/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
virtual

Set initial conditions.

Reimplemented from oomph::Problem.

◆ set_initial_condition() [10/10]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::set_initial_condition ( )
virtual

Set initial conditions.

Reimplemented from oomph::Problem.

◆ unsteady_run() [1/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

Perform run up to specified time t_max with given timestep dt.

416 {
417 
418  // Set value of epsilon
420 
421  // Deform the mesh/free surface
423 
424  // Initialise DocInfo object
425  DocInfo doc_info;
426 
427  // Set output directory
428  doc_info.set_directory("RESLT");
429 
430  // Initialise counter for solutions
431  doc_info.number()=0;
432 
433  // Open trace file
434  char filename[100];
435  sprintf(filename,"%s/trace.dat",doc_info.directory().c_str());
436  Trace_file.open(filename);
437 
438  // Initialise trace file
439  Trace_file << "time" << ", "
440  << "edge spine height" << ", "
441  << "contact angle left" << ", "
442  << "contact angle right" << ", " << std::endl;
443 
444  // Initialise timestep
445  initialise_dt(dt);
446 
447  // Set initial conditions
449 
450  // Determine number of timesteps
451  const unsigned n_timestep = unsigned(t_max/dt);
452 
453  // Open pvd file -- a wrapper for all the different
454  // vtu output files plus information about continuous time
455  // to facilitate animations in paraview
456  sprintf(filename,"%s/soln.pvd",doc_info.directory().c_str());
459 
460  // Doc initial solution
461  doc_solution(doc_info);
462 
463  // Increment counter for solutions
464  doc_info.number()++;
465 
466  // Timestepping loop
467  for(unsigned t=1;t<=n_timestep;t++)
468  {
469  // Output current timestep to screen
470  cout << "\nTimestep " << t << " of " << n_timestep << std::endl;
471 
472  // Take one fixed timestep
474 
475  // Doc solution
476  doc_solution(doc_info);
477 
478  // Increment counter for solutions
479  doc_info.number()++;
480  } // End of timestepping loop
481 
482  // write footer and close pvd file
485 
486 } // End of unsteady_run
void set_initial_condition()
Definition: rayleigh_instability.cc:122
void deform_free_surface(const double &epsilon)
Deform the mesh/free surface to a prescribed function.
Definition: rayleigh_instability.cc:160
void doc_solution(DocInfo &doc_info)
Doc the solution.
Definition: rayleigh_instability.cc:354
Definition: oomph_utilities.h:499
void set_directory(const std::string &directory)
Definition: oomph_utilities.cc:298
void initialise_dt(const double &dt)
Definition: problem.cc:13231
void unsteady_newton_solve(const double &dt)
Definition: problem.cc:10953
double Epsilon
Free surface cosine deform parameter.
Definition: rayleigh_instability.cc:71
void write_pvd_footer(std::ofstream &pvd_file)
Write the pvd file footer.
Definition: mesh.cc:9639
void write_pvd_header(std::ofstream &pvd_file)
Write the pvd file header.
Definition: mesh.cc:9615

References oomph::DocInfo::directory(), Global_Physical_Variables::Epsilon, oomph::SarahBL::epsilon, MergeRestartFiles::filename, oomph::DocInfo::number(), Global_Physical_Variables::Pvd_file, oomph::DocInfo::set_directory(), plotPSD::t, oomph::Problem_Parameter::Trace_file, oomph::ParaviewHelper::write_pvd_footer(), and oomph::ParaviewHelper::write_pvd_header().

◆ unsteady_run() [2/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [3/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [4/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [5/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [6/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [7/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [8/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [9/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [10/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const double t_max,
const double dt 
)

Do unsteady run up to maximum time t_max with given timestep dt.

◆ unsteady_run() [11/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const unsigned nstep)

Run an unsteady simulation with specified number of steps.

Unsteady run with specified number of steps.

485 {
486 
487  // Increase maximum residual
488  Problem::Max_residuals=500.0;
489 
490  //Distort the mesh
492  unsigned Nspine = Bulk_mesh_pt->nspine();
493  for(unsigned i=0;i<Nspine;i++)
494  {
495  double y_value = Bulk_mesh_pt->boundary_node_pt(0,i)->x(1);
496 
497  Bulk_mesh_pt->spine_pt(i)->height() =
500  }
501 
502  //Make sure to update
503  Bulk_mesh_pt->node_update();
504 
505  // Doc info object
506  DocInfo doc_info;
507 
508  // Set output directory
509  doc_info.set_directory("RESLT");
510  doc_info.number()=0;
511 
512  // Open trace file
513  char filename[100];
514  sprintf(filename,"%s/trace.dat",doc_info.directory().c_str());
515  Trace_file.open(filename);
516 
517  Trace_file << "VARIABLES=\"time\",";
518  Trace_file << "\"h<sub>left</sub>\",\"h<sub>right</sub>\"";
519  Trace_file << "\n";
520 
521  //Set value of dt
522  double dt = 0.1;
523 
524  //Initialise all time values
526 
527  //Doc initial solution
528  doc_solution(doc_info);
529 
530 //Loop over the timesteps
531  for(unsigned t=1;t<=nstep;t++)
532  {
533  cout << "TIMESTEP " << t << std::endl;
534 
535  //Take one fixed timestep
537 
538  // Doc solution
539  doc_info.number()++;
540  doc_solution(doc_info);
541  }
542 
543 }

References Global_Physical_Variables::Alpha, cos(), oomph::DocInfo::directory(), Global_Physical_Variables::Epsilon, oomph::SarahBL::epsilon, MergeRestartFiles::filename, Global_Physical_Variables::Film_Thickness, i, oomph::DocInfo::number(), oomph::DocInfo::set_directory(), plotPSD::t, and oomph::Problem_Parameter::Trace_file.

◆ unsteady_run() [12/12]

template<class ELEMENT , class TIMESTEPPER >
void InterfaceProblem< ELEMENT, TIMESTEPPER >::unsteady_run ( const unsigned nstep)

Run an unsteady simulation with specified number of steps.

Member Data Documentation

◆ Bulk_mesh_pt [1/8]

template<class ELEMENT , class TIMESTEPPER >
TwoLayerSpineMesh< ELEMENT > * InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt

Access function for the specific mesh.

Pointer to the specific bulk mesh.

The bulk fluid mesh, complete with spines and update information.

Pointer to the (specific) "bulk" mesh.

Referenced by InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem().

◆ Bulk_mesh_pt [2/8]

template<class ELEMENT , class TIMESTEPPER >
SingleLayerSpineMesh<ELEMENT>* InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt

Access function for the specific mesh.

The bulk fluid mesh, complete with spines and update information.

◆ Bulk_mesh_pt [3/8]

template<class ELEMENT , class TIMESTEPPER >
ElasticRefineableTwoLayerMesh<ELEMENT>* InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt
private

Pointer to the (specific) "bulk" mesh.

◆ Bulk_mesh_pt [4/8]

template<class ELEMENT , class TIMESTEPPER >
TwoLayerSpineMesh<ELEMENT>* InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt
private

Pointer to the specific bulk mesh.

◆ Bulk_mesh_pt [5/8]

template<class ELEMENT , class TIMESTEPPER >
AnnularSpineMesh<ELEMENT>* InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt
private

Pointer to bulk mesh.

◆ Bulk_mesh_pt [6/8]

template<class ELEMENT , class TIMESTEPPER >
MyHorizontalSingleLayerSpineMesh<ELEMENT>* InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt

Access function for the specific mesh.

◆ Bulk_mesh_pt [7/8]

template<class ELEMENT , class TIMESTEPPER >
ElasticRectangularQuadMesh<ELEMENT>* InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt

Access function for the specific mesh.

Pointer to the (specific) "bulk" mesh.

◆ Bulk_mesh_pt [8/8]

template<class ELEMENT , class TIMESTEPPER >
SingleLayerCubicSpineMesh<ELEMENT>* InterfaceProblem< ELEMENT, TIMESTEPPER >::Bulk_mesh_pt
private

Pointer to bulk mesh.

◆ Constitutive_law_pt

template<class ELEMENT , class TIMESTEPPER >
ConstitutiveLaw * InterfaceProblem< ELEMENT, TIMESTEPPER >::Constitutive_law_pt
private

Pointer to the constitutive law.

◆ Document_node_pt

template<class ELEMENT , class TIMESTEPPER >
Node * InterfaceProblem< ELEMENT, TIMESTEPPER >::Document_node_pt
private

Pointer to a node for documentation purposes.

Node for documentatin.

Referenced by InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem().

◆ Height

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::Height
private

Height of the domain.

◆ Interface_mesh_pt

template<class ELEMENT , class TIMESTEPPER >
Mesh * InterfaceProblem< ELEMENT, TIMESTEPPER >::Interface_mesh_pt

Mesh for the free surface (interface) elements.

Referenced by InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem().

◆ L_y

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::L_y
private

◆ Lr

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::Lr
private

Width of domain.

◆ Lx

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::Lx
private

Width of domain.

Axial lengths of domain.

◆ Ly

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::Ly
private

◆ R_max

template<class ELEMENT , class TIMESTEPPER >
double InterfaceProblem< ELEMENT, TIMESTEPPER >::R_max
private

Axial lengths of domain.

◆ Surface_mesh_pt

template<class ELEMENT , class TIMESTEPPER >
Mesh * InterfaceProblem< ELEMENT, TIMESTEPPER >::Surface_mesh_pt

Mesh for the interface elements.

The mesh that contains the free surface elements.

Pointer to the surface mes.

Pointer to the surface mesh.

Pointer to the "surface" mesh.

Referenced by InterfaceProblem< ELEMENT, TIMESTEPPER >::InterfaceProblem().

◆ Trace_file

template<class ELEMENT , class TIMESTEPPER >
ofstream InterfaceProblem< ELEMENT, TIMESTEPPER >::Trace_file
private

Trace file.


The documentation for this class was generated from the following files: