Public Member Functions | Private Attributes | List of all members
ElasticRingProblem< ELEMENT > Class Template Reference

Ring problem. More...

+ Inheritance diagram for ElasticRingProblem< ELEMENT >:

Public Member Functions

 ElasticRingProblem (const unsigned &n_element, bool &displ_control, bool &load_data_already_exists)
 Constructor for elastic ring problem. More...
 
OneDLagrangianMesh< ELEMENT > * mesh_pt ()
 Access function for the specific mesh. More...
 
void actions_after_newton_solve ()
 Update function is empty. More...
 
void actions_before_newton_solve ()
 Update function is empty. More...
 
void doc_solution (DocInfo &doc_info, ofstream &trace_file)
 Doc solution. More...
 
void parameter_study (DocInfo &doc_info)
 Perform the parameter study. More...
 
 ElasticRingProblem (const unsigned &n_element)
 Constructor: Pass umber of elements. More...
 
OneDLagrangianMesh< ELEMENT > * mesh_pt ()
 Access function for the specific mesh. More...
 
void actions_after_newton_solve ()
 Update function is empty. More...
 
void actions_before_newton_solve ()
 Update function is empty. More...
 
void doc_solution (DocInfo &doc_info, ofstream &trace_file)
 Doc solution. More...
 
void parameter_study (DocInfo &doc_info)
 Perform the parameter study. More...
 
 ElasticRingProblem (const unsigned &N, const double &L)
 Constructor for elastic ring problem. More...
 
OneDLagrangianMesh< ELEMENT > * mesh_pt ()
 Access function for the mesh. More...
 
void actions_after_newton_solve ()
 Update function is empty. More...
 
void actions_before_newton_solve ()
 Update function is empty. More...
 
void doc_solution (DocInfo &doc_info)
 Doc solution. More...
 
void unsteady_run ()
 Do unsteady run. More...
 
 ElasticRingProblem (const unsigned &n_element)
 Constructor: Number of elements. More...
 
OneDLagrangianMesh< ELEMENT > * mesh_pt ()
 Access function for the specific mesh. More...
 
void actions_after_newton_solve ()
 Update function is empty. More...
 
void actions_before_newton_solve ()
 Update function is empty. More...
 
void set_initial_conditions ()
 Setup initial conditions. More...
 
void doc_solution (DocInfo &doc_info)
 Doc solution. More...
 
void unsteady_run ()
 Do unsteady run. More...
 
void dump_it (ofstream &dump_file)
 
void restart (ifstream &restart_file)
 
- 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
 
virtual void actions_before_adapt ()
 
virtual void actions_after_adapt ()
 Actions that are to be performed after a mesh adaptation. More...
 
OomphCommunicator * communicator_pt ()
 access function to the oomph-lib communicator More...
 
const OomphCommunicator * communicator_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 ()
 
LinearSolver * mass_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...
 
Time * time_pt () const
 Return a pointer to the global time object (const version). More...
 
double & time ()
 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 TimeStepper * time_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...
 
double & minimum_dt ()
 Access function to min timestep in adaptive timestepping. More...
 
double & maximum_dt ()
 Access function to max timestep in adaptive timestepping. More...
 
unsigned & max_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...
 
double & max_residuals ()
 
bool & time_adaptive_newton_crash_on_solve_fail ()
 Access function for Time_adaptive_newton_crash_on_solve_fail. More...
 
double & newton_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...
 
double * global_dof_pt (const unsigned &i)
 
double & dof (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...
 
double * dof_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...
 
bool & use_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 Problem * make_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 ()
 
double * bifurcation_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 ()
 
int & sign_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 ()
 

Private Attributes

bool Displ_control
 Use displacement control? More...
 
GeomObject * Undef_geom_pt
 Pointer to geometric object that represents the undeformed shape. More...
 
unsigned Nbeam_element
 Number of elements in the beam mesh. More...
 
double Length
 Length of domain (in terms of the Lagrangian coordinates) More...
 
ELEMENT * Displ_control_elem_pt
 
Vector< double > S_displ_control
 At what local coordinate are we applying displacement control? More...
 
SolidInitialCondition * IC_pt
 Pointer to object that specifies the initial condition. More...
 
ofstream Trace_file
 Trace file for recording control data. More...
 
unsigned Validation_run_flag
 Flag for validation run: Default: 0 = no validation run. More...
 
bool Restart_flag
 Restart flag specified via command line? 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...
 
- Public Attributes inherited from oomph::Problem
bool Shut_up_in_newton_solve
 
- 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_convergence_check ()
 
virtual void actions_before_newton_step ()
 
virtual void actions_after_newton_step ()
 
virtual void actions_before_implicit_timestep ()
 
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)
 
double & dof_derivative (const unsigned &i)
 
double & dof_current (const unsigned &i)
 
virtual void set_initial_condition ()
 
virtual double global_temporal_error_norm ()
 
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 *, bool > Calculate_dparameter_analytic
 
bool Calculate_hessian_products_analytic
 
LinearAlgebraDistribution * Dof_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< double > Max_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< double > Dof_derivative
 Storage for the derivative of the problem variables wrt arc-length. More...
 
Vector< double > Dof_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
 
OomphCommunicator * Communicator_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 ElasticRingProblem< ELEMENT >

Ring problem.

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

////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////// Oscillating ring problem: Compare small-amplitude oscillations against analytical solution of the linearised equations.

Constructor & Destructor Documentation

◆ ElasticRingProblem() [1/4]

template<class ELEMENT >
ElasticRingProblem< ELEMENT >::ElasticRingProblem ( const unsigned &  n_element,
bool &  displ_control,
bool &  load_data_already_exists 
)

Constructor for elastic ring problem.

Constructor: Number of elements and flags for displ control and displacement control with existing data respectively.

150  :
151  Displ_control(displ_control),Nbeam_element(n_element)
152 {
153 
154  // Undeformed beam is an elliptical ring
155  Undef_geom_pt=new Ellipse(1.0,1.0);
156 
157  // Length of the doamin (in terms of the Lagrangian coordinates)
158  double length=2.0*atan(1.0);
159 
160  //Now create the (Lagrangian!) mesh
161  Problem::mesh_pt() =
162  new OneDLagrangianMesh<ELEMENT>(n_element,length,Undef_geom_pt);
163 
164  // Boundary condition:
165 
166  // Bottom:
167  unsigned ibound=0;
168  // No vertical displacement
169  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1);
170  // Infinite slope: Pin type 1 (slope) dof for displacement direction 0
171  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1,0);
172 
173  // Top:
174  ibound=1;
175  // No horizontal displacement
176  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(0);
177  // Zero slope: Pin type 1 (slope) dof for displacement direction 1
178  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1,1);
179 
180 
181  // Normal load incrementation
182  //===========================
183  if (!Displ_control)
184  {
185  // Create Data object whose one-and-only value contains the
186  // (in principle) adjustable load
187  Global_Physical_Variables::Pext_data_pt=new Data(1);
188 
189  //Pin the external pressure because it isn't actually adjustable.
190  Global_Physical_Variables::Pext_data_pt->pin(0);
191  }
192  // Displacement control
193  //=====================
194  else
195  {
196  // Choose element in which displacement control is applied: the last one
197  SolidFiniteElement* controlled_element_pt=
198  dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(Nbeam_element-1));
199 
200  // Fix the displacement in the vertical (1) direction...
201  unsigned controlled_direction=1;
202 
203  //... at right end of the control element
204  Vector<double> s_displ_control(1);
205  s_displ_control[0]=1.0;
206 
207  // Pointer to displacement control element
208  DisplacementControlElement* displ_control_el_pt;
209 
210  // Displacement control without previously existing load Data
211  //-----------------------------------------------------------
212  if (!load_data_already_exists)
213  {
214  // Build displacement control element
215  displ_control_el_pt=
216  new DisplacementControlElement(controlled_element_pt,
217  s_displ_control,
218  controlled_direction,
219  &Global_Physical_Variables::Xprescr);
220 
221  // The constructor of the DisplacementControlElement has created
222  // a new Data object whose one-and-only value contains the
223  // adjustable load: Use this Data object in the load function:
224  Global_Physical_Variables::Pext_data_pt=displ_control_el_pt->
225  displacement_control_load_pt();
226  }
227  // Demonstrate use of displacement control with some existing data
228  //----------------------------------------------------------------
229  else
230  {
231  // Create Data object whose one-and-only value contains the
232  // adjustable load
233  Global_Physical_Variables::Pext_data_pt=new Data(1);
234 
235  // Currently, nobody's "in charge of" this Data so it won't
236  // get included in any equation numbering schemes etc.
237  // --> declare it to be "global Data" for the Problem
238  // so the Problem is in charge and will perform such tasks.
239  add_global_data(Global_Physical_Variables::Pext_data_pt);
240 
241  // Build displacement control element and pass pointer to the
242  // already existing adjustable load Data.
243  displ_control_el_pt=
244  new DisplacementControlElement(controlled_element_pt,
245  s_displ_control,
246  controlled_direction,
247  &Global_Physical_Variables::Xprescr,
248  Global_Physical_Variables::Pext_data_pt);
249  }
250 
251  // Add the displacement-control element to the mesh
252  mesh_pt()->add_element_pt(displ_control_el_pt);
253  }
254 
255 
256 
257 
258  //Loop over the elements to set physical parameters etc.
259  for(unsigned i=0;i<Nbeam_element;i++)
260  {
261  //Cast to proper element type
262  ELEMENT *elem_pt = dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(i));
263 
264  // Set wall thickness
265  elem_pt->h_pt() = &Global_Physical_Variables::H;
266 
267  // Function that specifies load Vector
268  elem_pt->load_vector_fct_pt() = &Global_Physical_Variables::press_load;
269 
270  //Assign the undeformed beam shape
271  elem_pt->undeformed_beam_pt() = Undef_geom_pt;
272 
273  // Displacement control? If so, the load on *all* elements
274  // is affected by an unknown -- the external pressure, stored
275  // as the one-and-only value in a Data object: Add it to the
276  // elements' external Data.
277  if (Displ_control)
278  {
279  //The external pressure is external data for all elements
280  elem_pt->add_external_data(Global_Physical_Variables::Pext_data_pt);
281  }
282  }
283 
284  // Do equation numbering
285  cout << "# of dofs " << assign_eqn_numbers() << std::endl;
286 
287 } // end of constructor
i
int i
Definition: BiCGSTAB_step_by_step.cpp:9
ElasticRingProblem::Undef_geom_pt
GeomObject * Undef_geom_pt
Pointer to geometric object that represents the undeformed shape.
Definition: steady_ring.cc:132
ElasticRingProblem::mesh_pt
OneDLagrangianMesh< ELEMENT > * mesh_pt()
Access function for the specific mesh.
Definition: steady_ring.cc:109
ElasticRingProblem::Displ_control
bool Displ_control
Use displacement control?
Definition: steady_ring.cc:129
ElasticRingProblem::Nbeam_element
unsigned Nbeam_element
Number of elements in the beam mesh.
Definition: steady_ring.cc:135
oomph::Data
Definition: nodes.h:86
oomph::Data::pin
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:385
oomph::DisplacementControlElement
Definition: displacement_control_element.h:106
oomph::Ellipse
Definition: geom_objects.h:644
oomph::OneDLagrangianMesh
Definition: one_d_lagrangian_mesh.template.h:46
oomph::Problem::add_global_data
void add_global_data(Data *const &global_data_pt)
Definition: problem.h:1654
oomph::Problem::assign_eqn_numbers
unsigned long assign_eqn_numbers(const bool &assign_local_eqn_numbers=true)
Definition: problem.cc:1989
oomph::SolidFiniteElement
Definition: elements.h:3561
Eigen::bfloat16_impl::atan
EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bfloat16 atan(const bfloat16 &a)
Definition: BFloat16.h:636
Global_Physical_Variables::Xprescr
double Xprescr
Prescribed position (only used for displacement control)
Definition: steady_ring.cc:53
Global_Physical_Variables::press_load
void press_load(const Vector< double > &xi, const Vector< double > &x, const Vector< double > &N, Vector< double > &load)
Load function: Perturbation pressure to force non-axisymmetric deformation.
Definition: steady_ring.cc:64
Global_Physical_Variables::Pext_data_pt
Data * Pext_data_pt
Pointer to pressure load.
Definition: steady_ring.cc:60
Global_Physical_Variables::H
double H
Nondim thickness.
Definition: steady_ring.cc:50

References oomph::Problem::add_global_data(), oomph::Problem::assign_eqn_numbers(), Eigen::bfloat16_impl::atan(), ElasticRingProblem< ELEMENT >::Displ_control, Global_Physical_Variables::H, i, ElasticRingProblem< ELEMENT >::mesh_pt(), ElasticRingProblem< ELEMENT >::Nbeam_element, Global_Physical_Variables::Pext_data_pt, oomph::Data::pin(), Global_Physical_Variables::press_load(), ElasticRingProblem< ELEMENT >::Undef_geom_pt, and Global_Physical_Variables::Xprescr.

◆ ElasticRingProblem() [2/4]

template<class ELEMENT , class TIMESTEPPER >
ElasticRingProblem< ELEMENT, TIMESTEPPER >::ElasticRingProblem ( const unsigned &  n_element)

Constructor: Pass umber of elements.

Constructor for elastic ring problem.

152 {
153 
154  // Undeformed beam is an elliptical ring
155  Undef_geom_pt=new Ellipse(1.0,1.0);
156 
157  // Angle
158  double phi=MathematicalConstants::Pi/
159  double(Global_Physical_Variables::Nbuckl);
160 
161  // Length of the doamin (in terms of the Lagrangian coordinates)
162  double length=phi;
163 
164  // Actual number of elements
165  Nbeam_element=n_element;
166 
167  //Now create the (Lagrangian!) mesh
168  Problem::mesh_pt() =
169  new OneDLagrangianMesh<ELEMENT>(Nbeam_element,length,Undef_geom_pt);
170 
171  // "Bulk" beam element that the boundary condition element is attached to
172  int face_index=1;
173  ELEMENT* bulk_el_pt=dynamic_cast<ELEMENT*>(
174  mesh_pt()->element_pt(n_element-1));
175 
176  // Create boundary condition element
177  ClampedSlidingHermiteBeamBoundaryConditionElement*
178  bc_el_pt = new
179  ClampedSlidingHermiteBeamBoundaryConditionElement(
180  bulk_el_pt,face_index);
181 
182  // Set vectors that define the symmetry line
183  Vector<double> vector_to_symmetry_line(2);
184  vector_to_symmetry_line[0]=0.0;
185  vector_to_symmetry_line[1]=0.0;
186 
187  Vector<double> normal_to_symmetry_line(2);
188  normal_to_symmetry_line[0]=-sin(phi);
189  normal_to_symmetry_line[1]= cos(phi);
190 
191  bc_el_pt->set_symmetry_line(vector_to_symmetry_line,
192  normal_to_symmetry_line);
193 
194  // Add boundary condition element to mesh
195  mesh_pt()->add_element_pt(bc_el_pt);
196 
197 
198  // Boundary condition: Only at bottom
199  //-----------------------------------
200 
201  // Bottom:
202  unsigned ibound=0;
203 
204  // No vertical displacement
205  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1);
206 
207  // Infinite slope: Pin type 1 (slope) dof for displacement direction 0
208  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1,0);
209 
210 
211  // Displacement control
212  //---------------------
213 
214  // Choose element in which displacement control is applied: the first one
215  SolidFiniteElement* controlled_element_pt=
216  dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(0));
217 
218  // Fix the displacement in the horizontal direction...
219  unsigned controlled_direction=0;
220 
221  //... at left end of the control element
222  Vector<double> s_displ_control(1);
223  s_displ_control[0]=0.0;
224 
225  // Pointer to displacement control element
226  DisplacementControlElement* displ_control_el_pt;
227 
228 
229  // Build displacement control element
230  displ_control_el_pt=
231  new DisplacementControlElement(controlled_element_pt,
232  s_displ_control,
233  controlled_direction,
234  &Global_Physical_Variables::Xprescr);
235 
236  // The constructor of the DisplacementControlElement has created
237  // a new Data object whose one-and-only value contains the
238  // adjustable load: Use this Data object in the load function:
239  Global_Physical_Variables::Pext_data_pt=displ_control_el_pt->
240  displacement_control_load_pt();
241 
242 
243  // Add the displacement-control element to the mesh
244  mesh_pt()->add_element_pt(displ_control_el_pt);
245 
246  //Loop over the elements to set physical parameters etc.
247  for(unsigned i=0;i<Nbeam_element;i++)
248  {
249  //Cast to proper element type
250  ELEMENT *elem_pt = dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(i));
251 
252  // Set wall thickness
253  elem_pt->h_pt() = &Global_Physical_Variables::H;
254 
255  // Function that specifies load Vector
256  elem_pt->load_vector_fct_pt() = &Global_Physical_Variables::press_load;
257 
258  //Assign the undeformed beam shape
259  elem_pt->undeformed_beam_pt() = Undef_geom_pt;
260 
261  // Displacement control: The load on *all* elements
262  // is affected by an unknown -- the external pressure, stored
263  // as the one-and-only value in a Data object: Add it to the
264  // elements' external Data.
265  elem_pt->add_external_data(Global_Physical_Variables::Pext_data_pt);
266  }
267 
268 
269 
270  // Do equation numbering
271  cout << "# of dofs " << assign_eqn_numbers() << std::endl;
272 
273 } // end of constructor
cos
AnnoyingScalar cos(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:136
sin
AnnoyingScalar sin(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:137
oomph::ClampedSlidingHermiteBeamBoundaryConditionElement::set_symmetry_line
void set_symmetry_line(const Vector< double > &vector_to_symmetry_line, const Vector< double > &normal_to_symmetry_line)
Broken assignment operator.
Definition: beam_elements.h:548
BiharmonicTestFunctions2::Pi
double Pi
Definition: two_d_biharmonic.cc:235
Global_Physical_Variables::Nbuckl
unsigned Nbuckl
Buckling wavenumber.
Definition: steady_third_ring.cc:68

References cos(), Global_Physical_Variables::H, i, Global_Physical_Variables::Nbuckl, Global_Physical_Variables::Pext_data_pt, BiharmonicTestFunctions2::Pi, Global_Physical_Variables::press_load(), oomph::ClampedSlidingHermiteBeamBoundaryConditionElement::set_symmetry_line(), sin(), and Global_Physical_Variables::Xprescr.

◆ ElasticRingProblem() [3/4]

template<class ELEMENT , class TIMESTEPPER >
ElasticRingProblem< ELEMENT, TIMESTEPPER >::ElasticRingProblem ( const unsigned &  N,
const double &  L 
)

Constructor for elastic ring problem.

Constructor: Number of elements, length of domain, flag for setting Newmark IC directly or consistently

129  : Length(L)
130 {
131 
132  //Allocate the timestepper -- This constructs the time object as well
133  add_time_stepper_pt(new TIMESTEPPER());
134 
135  // Undeformed beam is an elliptical ring
136  Undef_geom_pt=new Ellipse(1.0,1.0);
137 
138  //Now create the (Lagrangian!) mesh
139  Problem::mesh_pt() = new OneDLagrangianMesh<ELEMENT>(
140  N,L,Undef_geom_pt,Problem::time_stepper_pt());
141 
142  // Boundary condition:
143 
144  // Bottom:
145  unsigned ibound=0;
146  // No vertical displacement
147  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1);
148  // Zero slope: Pin type 1 dof for displacement direction 0
149  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1,0);
150 
151  // Top:
152  ibound=1;
153  // No horizontal displacement
154  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(0);
155  // Zero slope: Pin type 1 dof for displacement direction 1
156  mesh_pt()->boundary_node_pt(ibound,0)->pin_position(1,1);
157 
158 
159  // Resize vector of local coordinates for control displacement
160  // (here only used to identify the point whose displacement we're
161  // tracing)
162  S_displ_control.resize(1);
163 
164  // Complete build of all elements so they are fully functional
165  // -----------------------------------------------------------
166 
167  // Find number of elements in mesh
168  unsigned Nelement = mesh_pt()->nelement();
169 
170  // Loop over the elements to set pointer to undeformed wall shape
171  for(unsigned i=0;i<Nelement;i++)
172  {
173  // Cast to proper element type
174  ELEMENT *elem_pt = dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(i));
175 
176  // Assign the undeformed surface
177  elem_pt->undeformed_beam_pt() = Undef_geom_pt;
178  }
179 
180  // Establish control displacment: (even though no displacement
181  // control is applied we still want to doc the displacement at the same point)
182 
183  // Choose element: (This is the last one)
184  Displ_control_elem_pt=dynamic_cast<ELEMENT*>(
185  mesh_pt()->element_pt(Nelement-1));
186 
187  // Fix/doc the displacement in the vertical (1) direction at right end of
188  // the control element
189  S_displ_control[0]=1.0;
190 
191  // Do equation numbering
192  cout << "# of dofs " << assign_eqn_numbers() << std::endl;
193 
194  // Geometric object that specifies the initial conditions
195  double eps_buckl=1.0e-2;
196  double HoR=dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(0))->h();
197  unsigned n_buckl=2;
198  unsigned imode=2;
199  GeomObject* ic_geom_object_pt=
200  new PseudoBucklingRing(eps_buckl,HoR,n_buckl,imode,
201  Problem::time_stepper_pt());
202 
203  // Setup object that specifies the initial conditions:
204  IC_pt = new SolidInitialCondition(ic_geom_object_pt);
205 
206 } // end of constructor
L
MatrixXd L
Definition: LLT_example.cpp:6
ElasticRingProblem::IC_pt
SolidInitialCondition * IC_pt
Pointer to object that specifies the initial condition.
Definition: lin_unsteady_ring.cc:114
ElasticRingProblem::S_displ_control
Vector< double > S_displ_control
At what local coordinate are we applying displacement control?
Definition: lin_unsteady_ring.cc:108
ElasticRingProblem::Displ_control_elem_pt
ELEMENT * Displ_control_elem_pt
Definition: lin_unsteady_ring.cc:105
ElasticRingProblem::Length
double Length
Length of domain (in terms of the Lagrangian coordinates)
Definition: lin_unsteady_ring.cc:101
oomph::GeomObject
Definition: geom_objects.h:101
oomph::Problem::add_time_stepper_pt
void add_time_stepper_pt(TimeStepper *const &time_stepper_pt)
Definition: problem.cc:1545
oomph::PseudoBucklingRing
Definition: pseudo_buckling_ring.h:56
oomph::SolidInitialCondition
Definition: elements.h:3496
N
@ N
Definition: constructor.cpp:22

References oomph::Problem::add_time_stepper_pt(), oomph::Problem::assign_eqn_numbers(), ElasticRingProblem< ELEMENT >::Displ_control_elem_pt, i, ElasticRingProblem< ELEMENT >::IC_pt, L, ElasticRingProblem< ELEMENT >::mesh_pt(), N, ElasticRingProblem< ELEMENT >::S_displ_control, and ElasticRingProblem< ELEMENT >::Undef_geom_pt.

◆ ElasticRingProblem() [4/4]

template<class ELEMENT >
ElasticRingProblem< ELEMENT >::ElasticRingProblem ( const unsigned &  n_element)

Constructor: Number of elements.

Member Function Documentation

◆ actions_after_newton_solve() [1/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_after_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

115 {}

◆ actions_after_newton_solve() [2/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_after_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

120 {}

◆ actions_after_newton_solve() [3/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_after_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

87 {}

◆ actions_after_newton_solve() [4/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_after_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

113 {}

◆ actions_before_newton_solve() [1/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_before_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

118 {}

◆ actions_before_newton_solve() [2/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_before_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

123 {}

◆ actions_before_newton_solve() [3/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_before_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

90 {}

◆ actions_before_newton_solve() [4/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::actions_before_newton_solve ( )
inlinevirtual

Update function is empty.

Reimplemented from oomph::Problem.

116 {}

◆ doc_solution() [1/4]

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

Doc solution.

Document solution.

215 {
216 
217  cout << "Doc-ing step " << doc_info.number()
218  << " for time " << time_stepper_pt()->time_pt()->time() << std::endl;
219 
220 
221  // Loop over all elements to get global kinetic and potential energy
222  unsigned Nelem=mesh_pt()->nelement();
223  double global_kin=0;
224  double global_pot=0;
225  double pot,kin;
226  for (unsigned ielem=0;ielem<Nelem;ielem++)
227  {
228  dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(ielem))->get_energy(pot,kin);
229  global_kin+=kin;
230  global_pot+=pot;
231  }
232 
233 
234  // Control displacement for initial condition object
235  Vector<double> xi_ctrl(1);
236  Vector<double> posn_ctrl(2);
237 
238  // Lagrangian coordinate of control point
239  xi_ctrl[0]=Displ_control_elem_pt->interpolated_xi(S_displ_control,0);
240 
241  // Get position
242  IC_pt->geom_object_pt()->position(xi_ctrl,posn_ctrl);
243 
244  // Write trace file: Time, control position, energies
245  Trace_file << time_pt()->time() << " "
246  << Displ_control_elem_pt->interpolated_x(S_displ_control,1)
247  << " " << global_pot << " " << global_kin
248  << " " << global_pot + global_kin
249  << " " << posn_ctrl[1]
250  << std::endl;
251 
252 
253  ofstream some_file;
254  char filename[100];
255 
256  // Number of plot points
257  unsigned npts=5;
258 
259  // Output solution
260  sprintf(filename,"%s/ring%i.dat",doc_info.directory().c_str(),
261  doc_info.number());
262  some_file.open(filename);
263  mesh_pt()->output(some_file,npts);
264  some_file.close();
265 
266  // Loop over all elements do dump out previous solutions
267  unsigned nsteps=time_stepper_pt()->nprev_values();
268  for (unsigned t=0;t<=nsteps;t++)
269  {
270  sprintf(filename,"%s/ring%i-%i.dat",doc_info.directory().c_str(),
271  doc_info.number(),t);
272  some_file.open(filename);
273  unsigned Nelem=mesh_pt()->nelement();
274  for (unsigned ielem=0;ielem<Nelem;ielem++)
275  {
276  dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(ielem))->
277  output(t,some_file,npts);
278  }
279  some_file.close();
280  }
281 
282  // Output for initial condition object
283  sprintf(filename,"%s/ic_ring%i.dat",doc_info.directory().c_str(),
284  doc_info.number());
285  some_file.open(filename);
286 
287  unsigned nplot=1+(npts-1)*mesh_pt()->nelement();
288  Vector<double> xi(1);
289  Vector<double> posn(2);
290  Vector<double> veloc(2);
291  Vector<double> accel(2);
292 
293  for (unsigned iplot=0;iplot<nplot;iplot++)
294  {
295  xi[0]=Length/double(nplot-1)*double(iplot);
296 
297  IC_pt->geom_object_pt()->position(xi,posn);
298  IC_pt->geom_object_pt()->dposition_dt(xi,1,veloc);
299  IC_pt->geom_object_pt()->dposition_dt(xi,2,accel);
300 
301  some_file << posn[0] << " " << posn[1] << " "
302  << xi[0] << " "
303  << veloc[0] << " " << veloc[1] << " "
304  << accel[0] << " " << accel[1] << " "
305  << sqrt(pow(posn[0],2)+pow(posn[1],2)) << " "
306  << sqrt(pow(veloc[0],2)+pow(veloc[1],2)) << " "
307  << sqrt(pow(accel[0],2)+pow(accel[1],2)) << " "
308  << std::endl;
309  }
310 
311  some_file.close();
312 } // end of doc solution
sqrt
AnnoyingScalar sqrt(const AnnoyingScalar &x)
Definition: AnnoyingScalar.h:134
ElasticRingProblem::Trace_file
ofstream Trace_file
Trace file for recording control data.
Definition: lin_unsteady_ring.cc:117
oomph::DocInfo::directory
std::string directory() const
Output directory.
Definition: oomph_utilities.h:524
oomph::DocInfo::number
unsigned & number()
Number used (e.g.) for labeling output files.
Definition: oomph_utilities.h:554
oomph::GeomObject::position
virtual void position(const Vector< double > &zeta, Vector< double > &r) const =0
Parametrised position on object at current time: r(zeta).
oomph::GeomObject::dposition_dt
virtual void dposition_dt(const Vector< double > &zeta, const unsigned &j, Vector< double > &drdt)
Definition: geom_objects.h:292
oomph::Problem::time_pt
Time *& time_pt()
Return a pointer to the global time object.
Definition: problem.h:1504
oomph::Problem::time_stepper_pt
TimeStepper *& time_stepper_pt()
Definition: problem.h:1524
oomph::SolidInitialCondition::geom_object_pt
GeomObject *& geom_object_pt()
Definition: elements.h:3511
oomph::TimeStepper::nprev_values
virtual unsigned nprev_values() const =0
Number of previous values available: 0 for static, 1 for BDF<1>,...
oomph::TimeStepper::time_pt
Time *const & time_pt() const
Access function for the pointer to time (const version)
Definition: timesteppers.h:572
oomph::Time::time
double & time()
Return the current value of the continuous time.
Definition: timesteppers.h:123
Eigen::bfloat16_impl::pow
EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bfloat16 pow(const bfloat16 &a, const bfloat16 &b)
Definition: BFloat16.h:625
MergeRestartFiles.filename
string filename
Definition: MergeRestartFiles.py:39
plotPSD.t
t
Definition: plotPSD.py:36
output
void output(std::ostream &outfile, const unsigned &nplot)
Overload output function.
Definition: overloaded_element_body.h:490

References oomph::DocInfo::directory(), MergeRestartFiles::filename, Global_Physical_Variables::Length, oomph::DocInfo::number(), output(), Eigen::ArrayBase< Derived >::pow(), sqrt(), plotPSD::t, and oomph::Problem_Parameter::Trace_file.

◆ doc_solution() [2/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::doc_solution ( DocInfo &  doc_info)

Doc solution.

◆ doc_solution() [3/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::doc_solution ( DocInfo &  doc_info,
ofstream &  trace_file 
)

Doc solution.

Document solution.

296 {
297  ofstream some_file;
298  char filename[100];
299 
300  // Number of plot points
301  unsigned npts=5;
302 
303  // Output solution
304  sprintf(filename,"%s/ring%i.dat",doc_info.directory().c_str(),
305  doc_info.number());
306  some_file.open(filename);
307  mesh_pt()->output(some_file,npts);
308  some_file.close();
309 
310  // Local coordinates of plot points at left and right end of domain
311  Vector<double> s_left(1);
312  s_left[0]=-1.0;
313  Vector<double> s_right(1);
314  s_right[0]=1.0;
315 
316  // Write trace file: Pressure, two radii
317  trace_file
318  << Global_Physical_Variables::Pext_data_pt->value(0)/
319  (pow(Global_Physical_Variables::H,3)/12.0)
320  << " "
321  << dynamic_cast<ELEMENT*>(mesh_pt()->
322  element_pt(0))->interpolated_x(s_left,0)
323  << " "
324  << dynamic_cast<ELEMENT*>
325  (mesh_pt()->element_pt(Nbeam_element-1))->interpolated_x(s_right,1)
326  << std::endl;
327 
328 
329 } // end of doc
oomph::Data::value
double value(const unsigned &i) const
Definition: nodes.h:293

References oomph::DocInfo::directory(), MergeRestartFiles::filename, Global_Physical_Variables::H, oomph::DocInfo::number(), Global_Physical_Variables::Pext_data_pt, Eigen::ArrayBase< Derived >::pow(), and oomph::Data::value().

◆ doc_solution() [4/4]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::doc_solution ( DocInfo &  doc_info,
ofstream &  trace_file 
)

Doc solution.

◆ dump_it()

template<class ELEMENT , class TIMESTEPPER >
void ElasticRingProblem< ELEMENT, TIMESTEPPER >::dump_it ( ofstream &  dump_file)

Dump problem-specific parameter values, then dump generic problem data.

327 {
328 
329  // Write Pcos
330  dump_file << Global_Physical_Variables::Pcos << " # Pcos" << std::endl;
331 
332  // Write validation run flag
333  dump_file << Validation_run_flag
334  << " # Validation run flag" << std::endl;
335 
336  // Call generic dump()
337  Problem::dump(dump_file);
338 
339 } // end of dump it
ElasticRingProblem::Validation_run_flag
unsigned Validation_run_flag
Flag for validation run: Default: 0 = no validation run.
Definition: unsteady_ring.cc:142
Global_Physical_Variables::Pcos
double Pcos
Perturbation pressure.
Definition: steady_ring.cc:56

References Global_Physical_Variables::Pcos.

◆ mesh_pt() [1/4]

template<class ELEMENT >
OneDLagrangianMesh<ELEMENT>* ElasticRingProblem< ELEMENT >::mesh_pt ( )
inline

Access function for the specific mesh.

110  {
111  return dynamic_cast<OneDLagrangianMesh<ELEMENT>*>(Problem::mesh_pt());
112  }

Referenced by ElasticRingProblem< ELEMENT >::ElasticRingProblem().

◆ mesh_pt() [2/4]

template<class ELEMENT >
OneDLagrangianMesh<ELEMENT>* ElasticRingProblem< ELEMENT >::mesh_pt ( )
inline

Access function for the specific mesh.

115  {
116  return dynamic_cast<OneDLagrangianMesh<ELEMENT>*>(Problem::mesh_pt());
117  }

◆ mesh_pt() [3/4]

template<class ELEMENT >
OneDLagrangianMesh<ELEMENT>* ElasticRingProblem< ELEMENT >::mesh_pt ( )
inline

Access function for the mesh.

82  {
83  return dynamic_cast<OneDLagrangianMesh<ELEMENT>*>(Problem::mesh_pt());
84  }

◆ mesh_pt() [4/4]

template<class ELEMENT >
OneDLagrangianMesh<ELEMENT>* ElasticRingProblem< ELEMENT >::mesh_pt ( )
inline

Access function for the specific mesh.

108  {
109  return dynamic_cast<OneDLagrangianMesh<ELEMENT>*>(Problem::mesh_pt());
110  }

◆ parameter_study() [1/2]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::parameter_study ( DocInfo &  doc_info)

Perform the parameter study.

Solver loop to perform parameter study.

Perturbation pressure

338 {
339 
340  //Open a trace file
341  char filename[100];
342  sprintf(filename,"%s/trace.dat",doc_info.directory().c_str());
343  ofstream trace_file(filename);
344  trace_file << "VARIABLES=\"p_e_x_t\",\"R_1\",\"R_2\"" << std::endl;
345  trace_file << "ZONE" << std::endl;
346 
347  //Output initial data
348  doc_solution(doc_info,trace_file);
349 
350  // Number of steps
351  unsigned nstep= 11; //51;
352 
353  // Increments
354  double displ_increment=1.0/double(nstep-1);
355  double p_buckl=3.00*pow(Global_Physical_Variables::H,3)/12.0;
356  double p_owc =5.22*pow(Global_Physical_Variables::H,3)/12.0;
357  double pext_increment=(p_owc-p_buckl)/double(nstep-1);
358 
359  // Set initial values for parameters that are to be incremented
360  Global_Physical_Variables::Xprescr=1.0+displ_increment;
361  Global_Physical_Variables::Pext_data_pt->set_value(0,p_buckl-pext_increment);
362 
363 
364 
365  // Without displacement control the Newton method converges very slowly
366  // as we approach the axisymmetric state: Allow more iterations before
367  // calling it a day...
368  if (Displ_control)
369  {
370  Max_newton_iterations=100;
371  }
372 
373 
374  // Downward loop over parameter incrementation with pcos>0
375  //--------------------------------------------------------
376 
378  Global_Physical_Variables::Pcos=1.0e-5;
379 
380 
381  // Downward loop over parameter incrementation
382  //---------------------------------------------
383  for(unsigned i=1;i<=nstep;i++)
384  {
385 
386  // Displacement control?
387  if (!Displ_control)
388  {
389  // Increment pressure
390  Global_Physical_Variables::external_pressure() += pext_increment;
391  }
392  else
393  {
394  // Increment control displacement
395  Global_Physical_Variables::Xprescr-=displ_increment;
396  }
397 
398  // Solve
399  newton_solve();
400 
401  // Doc solution
402  doc_info.number()++;
403  doc_solution(doc_info,trace_file);
404 
405  } // end of downward loop
406 
407  // Reset perturbation pressure
408  //----------------------------
409  Global_Physical_Variables::Pcos=0.0;
410 
411  // Set initial values for parameters that are to be incremented
412  Global_Physical_Variables::Xprescr-=displ_increment;
413  Global_Physical_Variables::external_pressure() += pext_increment;
414 
415  // Start new zone for tecplot
416  trace_file << "ZONE" << std::endl;
417 
418  // Upward loop over parameter incrementation
419  //------------------------------------------
420  for(unsigned i=nstep;i<2*nstep;i++)
421  {
422 
423  // Displacement control?
424  if (!Displ_control)
425  {
426  // Increment pressure
427  Global_Physical_Variables::external_pressure() -= pext_increment;
428  }
429  else
430  {
431  // Increment control displacement
432  Global_Physical_Variables::Xprescr+=displ_increment;
433  }
434 
435  // Solve
436  newton_solve();
437 
438  // Doc solution
439  doc_info.number()++;
440  doc_solution(doc_info,trace_file);
441 
442  }
443 
444 } // end of run
ElasticRingProblem::doc_solution
void doc_solution(DocInfo &doc_info, ofstream &trace_file)
Doc solution.
Definition: steady_ring.cc:294
oomph::Data::set_value
void set_value(const unsigned &i, const double &value_)
Definition: nodes.h:271
oomph::Problem::Max_newton_iterations
unsigned Max_newton_iterations
Maximum number of Newton iterations.
Definition: problem.h:599
oomph::Problem::newton_solve
void newton_solve()
Use Newton method to solve the problem.
Definition: problem.cc:8783
Global_Physical_Variables::external_pressure
double & external_pressure()
Access function to value of external pressure.
Definition: steady_ring.cc:79

References oomph::DocInfo::directory(), Global_Physical_Variables::external_pressure(), MergeRestartFiles::filename, Global_Physical_Variables::H, i, oomph::Locate_zeta_helpers::Max_newton_iterations, oomph::DocInfo::number(), Global_Physical_Variables::Pcos, Global_Physical_Variables::Pext_data_pt, Eigen::ArrayBase< Derived >::pow(), oomph::Data::set_value(), and Global_Physical_Variables::Xprescr.

◆ parameter_study() [2/2]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::parameter_study ( DocInfo &  doc_info)

Perform the parameter study.

◆ restart()

template<class ELEMENT , class TIMESTEPPER >
void ElasticRingProblem< ELEMENT, TIMESTEPPER >::restart ( ifstream &  restart_file)

Read problem-specific parameter values, then recover generic problem data.

348 {
349 
350  string input_string;
351 
352  // Read line up to termination sign
353  getline(restart_file,input_string,'#');
354  // Ignore rest of line
355  restart_file.ignore(80,'\n');
356  // Read in pcos
357  Global_Physical_Variables::Pcos=atof(input_string.c_str());
358 
359  // Read line up to termination sign
360  getline(restart_file,input_string,'#');
361  // Ignore rest of line
362  restart_file.ignore(80,'\n');
363  // Read in Long run flag
364  Validation_run_flag=
365  unsigned(atof(input_string.c_str()));
366 
367  // Call generic read()
368  Problem::read(restart_file);
369 
370 } // end of restart
Eigen::TensorSycl::internal::read
static EIGEN_DEVICE_FUNC EIGEN_STRONG_INLINE std::enable_if_t< PacketLoad, PacketType > read(const TensorMapper &tensorMapper, const StorageIndex &NCIndex, const StorageIndex &CIndex, const StorageIndex &ld)
read, a template function used for loading the data from global memory. This function is used to guar...
Definition: TensorContractionSycl.h:162

References Global_Physical_Variables::Pcos, and Eigen::TensorSycl::internal::read().

◆ set_initial_conditions()

template<class ELEMENT , class TIMESTEPPER >
void ElasticRingProblem< ELEMENT, TIMESTEPPER >::set_initial_conditions

Setup initial conditions.

Setup initial conditions – either restart from solution specified via command line or impulsive start.

381 {
382 
383  // No restart file --> impulsive start from initial configuration
384  // assigned in the Lagrangian mesh.
385  if (!Restart_flag)
386  {
387  // Set initial timestep
388  double dt=1.0;
389 
390  // Assign impulsive start
391  assign_initial_values_impulsive(dt);
392  }
393  // Restart file specified via command line
394  else
395  {
396  // Try to open restart file
397  ifstream* restart_file_pt=
398  new ifstream(CommandLineArgs::Argv[2],ios_base::in);
399  if (restart_file_pt!=0)
400  {
401  cout << "Have opened " << CommandLineArgs::Argv[2] <<
402  " for restart. " << std::endl;
403  }
404  else
405  {
406  std::ostringstream error_stream;
407  error_stream << "ERROR while trying to open "
408  << CommandLineArgs::Argv[2]
409  << " for restart." << std::endl;
410 
411  throw OomphLibError(error_stream.str(),
412  OOMPH_CURRENT_FUNCTION,
413  OOMPH_EXCEPTION_LOCATION);
414  }
415 
416  // Read restart data:
417  restart(*restart_file_pt);
418 
419  }
420 
421 } // end of set ic
ElasticRingProblem::restart
void restart(ifstream &restart_file)
Definition: unsteady_ring.cc:347
ElasticRingProblem::Restart_flag
bool Restart_flag
Restart flag specified via command line?
Definition: unsteady_ring.cc:145
oomph::OomphLibError
Definition: oomph_definitions.h:222
oomph::Problem::assign_initial_values_impulsive
void assign_initial_values_impulsive()
Definition: problem.cc:11499
oomph::CommandLineArgs::Argv
char ** Argv
Arguments themselves.
Definition: oomph_utilities.cc:410
OOMPH_EXCEPTION_LOCATION
#define OOMPH_EXCEPTION_LOCATION
Definition: oomph_definitions.h:61
OOMPH_CURRENT_FUNCTION
#define OOMPH_CURRENT_FUNCTION
Definition: oomph_definitions.h:86

References oomph::CommandLineArgs::Argv, OOMPH_CURRENT_FUNCTION, and OOMPH_EXCEPTION_LOCATION.

◆ unsteady_run() [1/2]

template<class ELEMENT , class TIMESTEPPER >
void ElasticRingProblem< ELEMENT, TIMESTEPPER >::unsteady_run

Do unsteady run.

Solver loop to perform unsteady run.

Label for output

321 {
322 
324  DocInfo doc_info;
325 
326  // Output directory
327  doc_info.set_directory("RESLT");
328 
329  // Step number
330  doc_info.number()=0;
331 
332  // Set up trace file
333  char filename[100];
334  sprintf(filename,"%s/trace_ring.dat",doc_info.directory().c_str());
335  Trace_file.open(filename);
336 
337  Trace_file << "VARIABLES=\"time\",\"R<sub>ctrl</sub>\",\"E<sub>pot</sub>\"";
338  Trace_file << ",\"E<sub>kin</sub>\",\"E<sub>kin</sub>+E<sub>pot</sub>\"";
339  Trace_file << ",\"<sub>exact</sub>R<sub>ctrl</sub>\""
340  << std::endl;
341 
342  // Number of steps
343  unsigned nstep=600;
344  if (Global_Physical_Variables::Long_run_flag==0) {nstep=10;}
345 
346  // Initial timestep
347  double dt=1.0;
348 
349  // Ratio for timestep reduction
350  double timestep_ratio=1.0;
351  if (Global_Physical_Variables::Fixed_timestep_flag==0) {timestep_ratio=0.995;}
352 
353  // Number of previous timesteps stored
354  unsigned ndt=time_stepper_pt()->time_pt()->ndt();
355 
356  // Setup vector of "previous" timesteps
357  Vector<double> dt_prev(ndt);
358  dt_prev[0]=dt;
359  for (unsigned i=1;i<ndt;i++)
360  {
361  dt_prev[i]=dt_prev[i-1]/timestep_ratio;
362  }
363 
364  // Initialise the history of previous timesteps
365  time_pt()->initialise_dt(dt_prev);
366 
367  // Initialise time
368  double time0=10.0;
369  time_pt()->time()=time0;
370 
371  // Setup analytical initial condition?
372  if (Global_Physical_Variables::Consistent_newmark_ic)
373  {
374  // Note: Time has been scaled on intrinsic timescale so
375  // we don't need to specify a multiplier for the inertia
376  // terms (the default assignment of 1.0 is OK)
377  SolidMesh::Solid_IC_problem.
378  set_newmark_initial_condition_consistently(
379  this,mesh_pt(),static_cast<TIMESTEPPER*>(time_stepper_pt()),IC_pt,dt);
380  }
381  else
382  {
383  SolidMesh::Solid_IC_problem.
384  set_newmark_initial_condition_directly(
385  this,mesh_pt(),static_cast<TIMESTEPPER*>(time_stepper_pt()),IC_pt,dt);
386  }
387 
388  //Output initial data
389  doc_solution(doc_info);
390 
391  // Time integration loop
392  for(unsigned i=1;i<=nstep;i++)
393  {
394  // Solve
395  unsteady_newton_solve(dt);
396 
397  // Doc solution
398  doc_info.number()++;
399  doc_solution(doc_info);
400 
401  // Reduce timestep
402  if (time_pt()->time()<100.0) {dt=timestep_ratio*dt;}
403  }
404 
405 } // end of unsteady run
oomph::DocInfo
Definition: oomph_utilities.h:499
oomph::DocInfo::set_directory
void set_directory(const std::string &directory)
Definition: oomph_utilities.cc:298
oomph::Problem::unsteady_newton_solve
void unsteady_newton_solve(const double &dt)
Definition: problem.cc:10953
oomph::Problem::time
double & time()
Return the current value of continuous time.
Definition: problem.cc:11531
oomph::Time::initialise_dt
void initialise_dt(const double &dt_)
Set all timesteps to the same value, dt.
Definition: timesteppers.h:99
oomph::Time::ndt
unsigned ndt() const
Return the number of timesteps stored.
Definition: timesteppers.h:129
Global_Physical_Variables::Fixed_timestep_flag
unsigned Fixed_timestep_flag
Flag for fixed timestep: Default = fixed timestep.
Definition: lin_unsteady_ring.cc:52
Global_Physical_Variables::Long_run_flag
unsigned Long_run_flag
Flag for long/short run: Default = perform long run.
Definition: lin_unsteady_ring.cc:49
Global_Physical_Variables::Consistent_newmark_ic
bool Consistent_newmark_ic
Definition: lin_unsteady_ring.cc:56

References Global_Physical_Variables::Consistent_newmark_ic, oomph::DocInfo::directory(), MergeRestartFiles::filename, Global_Physical_Variables::Fixed_timestep_flag, i, Global_Physical_Variables::Long_run_flag, oomph::DocInfo::number(), oomph::DocInfo::set_directory(), and oomph::Problem_Parameter::Trace_file.

◆ unsteady_run() [2/2]

template<class ELEMENT >
void ElasticRingProblem< ELEMENT >::unsteady_run ( )

Do unsteady run.

Member Data Documentation

◆ Displ_control

template<class ELEMENT >
bool ElasticRingProblem< ELEMENT >::Displ_control
private

Use displacement control?

Referenced by ElasticRingProblem< ELEMENT >::ElasticRingProblem().

◆ Displ_control_elem_pt

template<class ELEMENT >
ELEMENT* ElasticRingProblem< ELEMENT >::Displ_control_elem_pt
private

In which element are we applying displacement control? (here only used for doc of radius)

Referenced by ElasticRingProblem< ELEMENT >::ElasticRingProblem().

◆ IC_pt

template<class ELEMENT >
SolidInitialCondition* ElasticRingProblem< ELEMENT >::IC_pt
private

Pointer to object that specifies the initial condition.

Referenced by ElasticRingProblem< ELEMENT >::ElasticRingProblem().

◆ Length

template<class ELEMENT >
double ElasticRingProblem< ELEMENT >::Length
private

Length of domain (in terms of the Lagrangian coordinates)

◆ Nbeam_element

template<class ELEMENT >
unsigned ElasticRingProblem< ELEMENT >::Nbeam_element
private

Number of elements in the beam mesh.

Referenced by ElasticRingProblem< ELEMENT >::ElasticRingProblem().

◆ Restart_flag

template<class ELEMENT >
bool ElasticRingProblem< ELEMENT >::Restart_flag
private

Restart flag specified via command line?

◆ S_displ_control

template<class ELEMENT >
Vector<double> ElasticRingProblem< ELEMENT >::S_displ_control
private

At what local coordinate are we applying displacement control?

Referenced by ElasticRingProblem< ELEMENT >::ElasticRingProblem().

◆ Trace_file

template<class ELEMENT >
ofstream ElasticRingProblem< ELEMENT >::Trace_file
private

Trace file for recording control data.

◆ Undef_geom_pt

template<class ELEMENT >
GeomObject * ElasticRingProblem< ELEMENT >::Undef_geom_pt
private

Pointer to geometric object that represents the undeformed shape.

Referenced by ElasticRingProblem< ELEMENT >::ElasticRingProblem().

◆ Validation_run_flag

template<class ELEMENT >
unsigned ElasticRingProblem< ELEMENT >::Validation_run_flag
private

Flag for validation run: Default: 0 = no validation run.

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