ShellProblem< ELEMENT > Class Template Reference
+ Inheritance diagram for ShellProblem< ELEMENT >:

Public Member Functions

 ShellProblem (const unsigned &nx, const unsigned &ny, const double &lx, const double &ly)
 Constructor. More...
 
 ~ShellProblem ()
 Destructor: delete mesh, geometric object. More...
 
ShellMesh< ELEMENT > * solid_mesh_pt ()
 Overload Access function for the mesh. More...
 
void actions_after_newton_solve ()
 Actions after solve empty. More...
 
void actions_before_newton_solve ()
 Actions before solve empty. More...
 
void actions_after_distribute ()
 
void solve ()
 
 ShellProblem (const unsigned &nx, const unsigned &ny, const double &lx, const double &ly, const unsigned &problem_id)
 Constructor. More...
 
CircularCylindricalShellMesh< ELEMENT > * mesh_pt ()
 Overload Access function for the mesh. More...
 
void actions_after_newton_solve ()
 Actions after solve empty. More...
 
void actions_before_newton_solve ()
 Actions before solve empty. More...
 
void create_bc_elements (const unsigned &b)
 Create clamping bc face elements on the b-th boundary of the Mesh. More...
 
void solve (string &dir_name)
 Do parameter study. More...
 
 ShellProblem (const unsigned &nx, const unsigned &ny, const double &lx, const double &ly)
 Constructor. More...
 
ShellMesh< ELEMENT > * mesh_pt ()
 Overload Access function for the mesh. More...
 
void actions_after_newton_solve ()
 Actions after solve empty. More...
 
void actions_before_newton_solve ()
 Actions before solve empty. More...
 
void solve ()
 
 ShellProblem (const unsigned &nx, const unsigned &ny, const double &lx, const double &ly)
 Constructor. More...
 
 ~ShellProblem ()
 Destructor: delete mesh, geometric object. More...
 
ShellMesh< ELEMENT > * mesh_pt ()
 Overload Access function for the mesh. More...
 
void actions_after_newton_solve ()
 Actions after solve empty. More...
 
void actions_before_newton_solve ()
 Actions before solve empty. More...
 
void solve ()
 
 ShellProblem (const unsigned &nx, const unsigned &ny, const double &lx, const double &ly)
 Constructor. More...
 
CircularCylindricalShellMesh< ELEMENT > * mesh_pt ()
 Overload Access function for the mesh. More...
 
void actions_after_newton_solve ()
 Actions after solve empty. More...
 
void actions_before_newton_solve ()
 Actions before solve empty. More...
 
void run_it ()
 Do parameter study. More...
 
void doc_solution (DocInfo &doc_info, ofstream &trace_file)
 Doc solution. 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
 
virtual void actions_before_adapt ()
 
virtual void actions_after_adapt ()
 Actions that are to be performed after a mesh adaptation. More...
 
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 ()
 

Private Attributes

GeomObjectUndeformed_midplane_pt
 Pointer to GeomObject that specifies the undeformed midplane. More...
 
NodeTrace_node_pt
 First trace node. More...
 
NodeTrace_node2_pt
 Second trace node. More...
 
MeshSolid_mesh_pt
 Pointer to the solid mesh. More...
 
MeshDisplacement_control_mesh_pt
 Pointer to the mesh that contains the displacement control element. More...
 
unsigned Nshell
 Number of shell elements. 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)
 
doubledof_derivative (const unsigned &i)
 
doubledof_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 *, 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...
 

Constructor & Destructor Documentation

◆ ShellProblem() [1/5]

template<class ELEMENT >
ShellProblem< ELEMENT >::ShellProblem ( const unsigned nx,
const unsigned ny,
const double lx,
const double ly 
)

Constructor.

410 {
411  Mesh::Suppress_warning_about_empty_mesh_level_time_stepper_function=true;
412 
414  //Create the undeformed midplane object
415  Undeformed_midplane_pt = new EllipticalTube(1.0,1.0);
416 
417  //Now create the mesh
419 
420  //Set the undeformed positions in the mesh
421  solid_mesh_pt()->assign_undeformed_positions(Undeformed_midplane_pt);
422 
423  //Reorder the elements, since I know what's best for them....
424  solid_mesh_pt()->element_reorder();
425 
426  //Apply boundary conditions to the ends of the tube
427  unsigned n_ends = solid_mesh_pt()->nboundary_node(1);
428  //Loop over the node
429  for(unsigned i=0;i<n_ends;i++)
430  {
431  //Pin in the axial direction (prevents rigid body motions)
432  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(2);
433  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(2);
434  //Derived conditions
435  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(2,2);
436  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(2,2);
437 
438  //------------------CLAMPING CONDITIONS----------------------
439  //------Pin positions in the transverse directions-----------
440  // Comment these out to get the ring case
441  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(0);
442  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(0);
443  //Derived conditions
444  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(2,0);
445  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(2,0);
446 
447  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(1);
448  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(1);
449  //Derived conditions
450  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(2,1);
451  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(2,1);
452  //----------------------------------------------------------
453 
454  // Set the axial gradients of the transverse coordinates to be
455  // zero --- need to be enforced for ring or tube buckling
456  //Pin dx/dz and dy/dz
457  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(1,0);
458  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(1,1);
459  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(1,0);
460  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(1,1);
461  //Derived conditions
462  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(3,0);
463  solid_mesh_pt()->boundary_node_pt(1,i)->pin_position(3,1);
464  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(3,0);
465  solid_mesh_pt()->boundary_node_pt(3,i)->pin_position(3,1);
466  }
467 
468  //Now loop over the sides and apply symmetry conditions
469  unsigned n_side = solid_mesh_pt()->nboundary_node(0);
470  for(unsigned i=0;i<n_side;i++)
471  {
472  //At the side where theta is 0, pin in the y direction
473  solid_mesh_pt()->boundary_node_pt(0,i)->pin_position(1);
474  //Derived condition
475  solid_mesh_pt()->boundary_node_pt(0,i)->pin_position(1,1);
476  //Pin dx/dtheta and dz/dtheta
477  solid_mesh_pt()->boundary_node_pt(0,i)->pin_position(2,0);
478  solid_mesh_pt()->boundary_node_pt(0,i)->pin_position(2,2);
479  //Pin the mixed derivative
480  solid_mesh_pt()->boundary_node_pt(0,i)->pin_position(3,0);
481  solid_mesh_pt()->boundary_node_pt(0,i)->pin_position(3,2);
482 
483  //At the side when theta is 0.5pi pin in the x direction
484  solid_mesh_pt()->boundary_node_pt(2,i)->pin_position(0);
485  //Derived condition
486  solid_mesh_pt()->boundary_node_pt(2,i)->pin_position(1,0);
487  //Pin dy/dtheta and dz/dtheta
488  solid_mesh_pt()->boundary_node_pt(2,i)->pin_position(2,1);
489  solid_mesh_pt()->boundary_node_pt(2,i)->pin_position(2,2);
490  //Pin the mixed derivative
491  solid_mesh_pt()->boundary_node_pt(2,i)->pin_position(3,1);
492  solid_mesh_pt()->boundary_node_pt(2,i)->pin_position(3,2);
493 
494  //Set an initial kick to make sure that we hop onto the
495  //non-axisymmetric branch
496  if((i>1) && (i<n_side-1))
497  {
498  solid_mesh_pt()->boundary_node_pt(0,i)->x(0) += 0.05;
499  solid_mesh_pt()->boundary_node_pt(2,i)->x(1) -= 0.1;
500  }
501  }
502 
503 
504  // Setup displacement control
505  //---------------------------
506 
507 
508 
509 // //Setup displacement control
510 // //Fix the displacement at the mid-point of the tube in the "vertical"
511 // //(y) direction.
512 // //Set the displacement control element (located halfway along the tube)
513 // Disp_ctl_element_pt = dynamic_cast<ELEMENT*>(solid_mesh_pt()->element_pt(3*Ny-1));
514 // //The midpoint of the tube is located exactly half-way along the element
515 // Vector<double> s(2); s[0] = 1.0; s[1] = 0.0; //s[1] = 0.5
516 // //Fix the displacement at this point in the y (1) direction
517 // Disp_ctl_element_pt->fix_displacement_for_displacement_control(s,1);
518 // //Set the pointer to the prescribed position
519 // Disp_ctl_element_pt->prescribed_position_pt() = &Prescribed_y;
520 
521 
522 
523  // Choose element in which displacement control is applied: This
524  // one is located halfway along the tube)
525  SolidFiniteElement* controlled_element_pt=
526  dynamic_cast<ELEMENT*>(solid_mesh_pt()->element_pt(3*ny-1));
527 
528  // Fix the displacement in the y (1) direction...
529  unsigned controlled_direction=1;
530 
531  // Local coordinate of the control point within the controlled element
532  Vector<double> s_displ_control(2);
533  s_displ_control[0]=1.0;
534  s_displ_control[1]=0.0;
535 
536  // Pointer to displacement control element
537  DisplacementControlElement* displ_control_el_pt;
538 
539  // Build displacement control element
540  displ_control_el_pt=
541  new DisplacementControlElement(controlled_element_pt,
542  s_displ_control,
543  controlled_direction,
545 
546  // The constructor of the DisplacementControlElement has created
547  // a new Data object whose one-and-only value contains the
548  // adjustable load: Use this Data object in the load function:
549  Global_Physical_Variables::Pext_data_pt=displ_control_el_pt->
550  displacement_control_load_pt();
551 
552  // Add the displacement-control element to the mesh
554  Displacement_control_mesh_pt->add_element_pt(displ_control_el_pt);
555 
556  //This mesh must be retained as halo on all processors
557 #ifdef OOMPH_HAS_MPI
558  Displacement_control_mesh_pt->set_keep_all_elements_as_halos();
559 #endif
560 
561 
562  // Complete build of shell elements
563  //---------------------------------
564 
565  //Find number of shell elements in mesh
566  unsigned n_element = nx*ny;
567 
568  //Explicit pointer to first element in the mesh
569  ELEMENT* first_el_pt = dynamic_cast<ELEMENT*>(solid_mesh_pt()->element_pt(0));
570 
571  //Loop over the elements
572  for(unsigned e=0;e<n_element;e++)
573  {
574  //Cast to a shell element
575  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(solid_mesh_pt()->element_pt(e));
576 
577  //Set the load function
578  el_pt->load_vector_fct_pt() = & Global_Physical_Variables::press_load;
579 
580  //Set the undeformed surface
581  el_pt->undeformed_midplane_pt() = Undeformed_midplane_pt;
582 
583  //The external pressure is external data for all elements
584  el_pt->add_external_data(Global_Physical_Variables::Pext_data_pt);
585 
586  //Pre-compute the second derivatives wrt Lagrangian coordinates
587  //for the first element only
588  if(e==0)
589  {
590  el_pt->pre_compute_d2shape_lagrangian_at_knots();
591  }
592  //Otherwise set the values to be the same as those in the first element
593  //this is OK because the Lagrangian mesh is uniform.
594  //
595  else
596  {
597  el_pt->set_dshape_lagrangian_stored_from_element(first_el_pt);
598  }
599  }
600 
601  //Set pointers to two trace nodes, used for output
602  //in the middle of the shell so available on both processors!
603  Trace_node_pt = solid_mesh_pt()->finite_element_pt(2*ny-1)->node_pt(3);
604  Trace_node2_pt = solid_mesh_pt()->finite_element_pt(ny)->node_pt(1);
605 
606  //Now add the two submeshes
609  this->build_global_mesh();
610 
611  // Do equation numbering
612  cout << std::endl;
613  cout << "------------------DISPLACEMENT CONTROL--------------------"
614  << std::endl;
615  cout << "# of dofs " << assign_eqn_numbers() << std::endl;
616  cout << "----------------------------------------------------------"
617  << std::endl;
618  cout << std::endl;
619 
620 }
int i
Definition: BiCGSTAB_step_by_step.cpp:9
Array< double, 1, 3 > e(1./3., 0.5, 2.)
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:89
Node * Trace_node2_pt
Second trace node.
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:392
GeomObject * Undeformed_midplane_pt
Pointer to GeomObject that specifies the undeformed midplane.
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:386
Node * Trace_node_pt
First trace node.
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:389
ShellMesh< ELEMENT > * solid_mesh_pt()
Overload Access function for the mesh.
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:334
Mesh * Solid_mesh_pt
Pointer to the solid mesh.
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:395
Mesh * Displacement_control_mesh_pt
Pointer to the mesh that contains the displacement control element.
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:398
Definition: displacement_control_element.h:106
Definition: geom_objects.h:1131
Definition: mesh.h:67
void add_element_pt(GeneralisedElement *const &element_pt)
Add a (pointer to) an element to the mesh.
Definition: mesh.h:617
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
bool Use_continuation_timestepper
Boolean to control original or new storage of dof stuff.
Definition: problem.h:758
Definition: elements.h:3561
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
double Prescribed_y
Prescribed position of control point.
Definition: mpi/distribution/clamped_shell/clamped_shell_with_arclength_cont.cc:54
Data * Pext_data_pt
Pointer to pressure load.
Definition: steady_ring.cc:60
const double ly
Definition: ConstraintElementsUnitTest.cpp:34
const double lx
Definition: ConstraintElementsUnitTest.cpp:33
const unsigned nx
Definition: ConstraintElementsUnitTest.cpp:30
const unsigned ny
Definition: ConstraintElementsUnitTest.cpp:31

References oomph::Mesh::add_element_pt(), e(), i, Mesh_Parameters::lx, Mesh_Parameters::ly, Mesh_Parameters::nx, Mesh_Parameters::ny, Global_Physical_Variables::Pext_data_pt, Global_Physical_Variables::Prescribed_y, and Global_Physical_Variables::press_load().

◆ ~ShellProblem() [1/2]

template<class ELEMENT >
ShellProblem< ELEMENT >::~ShellProblem ( )
inline

◆ ShellProblem() [2/5]

template<class ELEMENT >
ShellProblem< ELEMENT >::ShellProblem ( const unsigned nx,
const unsigned ny,
const double lx,
const double ly,
const unsigned problem_id 
)

Constructor.

138 {
139  //Create the undeformed midplane object
140  Undeformed_midplane_pt = new EllipticalTube(1.0,1.0);
141 
142  //Now create the mesh
143  Problem::mesh_pt() = new CircularCylindricalShellMesh<ELEMENT>(nx,ny,lx,ly);
144 
145  // Store number of genuine shell elements
146  Nshell=mesh_pt()->nelement();
147 
148  //Set the undeformed positions in the mesh
149  mesh_pt()->assign_undeformed_positions(Undeformed_midplane_pt);
150 
151  //Reorder the elements, since I know what's best for them....
152  mesh_pt()->element_reorder();
153 
154 
155 
156  //Loop over the nodes at the ends of the tube and pin the positions
157  //-----------------------------------------------------------------
158  // (All constraints apply at all nodes on boundaries 1 and 3)
159  //-----------------------------------------------------------
160  unsigned n_node_at_ends = mesh_pt()->nboundary_node(1);
161  for(unsigned j=0;j<n_node_at_ends;j++)
162  {
163 
164  // Loop over all three displacement components
165  for (unsigned i=0;i<3;i++)
166  {
167  // Pin positions for all s_1 (along tube circumference)
168  mesh_pt()->boundary_node_pt(1,j)->pin_position(i);
169  mesh_pt()->boundary_node_pt(3,j)->pin_position(i);
170 
171  // ...implying that dr/ds_1=0, too: [pin_position() takes type (2) and
172  // direction (i)]
173  mesh_pt()->boundary_node_pt(1,j)->pin_position(2,i);
174  mesh_pt()->boundary_node_pt(3,j)->pin_position(2,i);
175  }
176 
177 
178  // Use symmetry conditions at ends?
179  //--------------------------------
180  if (problem_id==2)
181  {
182  // Pin the derivative (into the shell) of the x-position for all
183  // s_1 (along tube circumference)[pin_position() takes type (1)
184  // and direction(0)]
185  mesh_pt()->boundary_node_pt(1,j)->pin_position(1,0);
186  mesh_pt()->boundary_node_pt(3,j)->pin_position(1,0);
187 
188  // ...implying that d^2x/ds_0 ds_1=0, too: [pin_position() takes
189  // type (3) and direction (0)]
190  mesh_pt()->boundary_node_pt(1,j)->pin_position(3,0);
191  mesh_pt()->boundary_node_pt(3,j)->pin_position(3,0);
192 
193  // Pin the derivative (into the shell) of the y-position for all
194  // s_1 (along tube circumference)[pin_position() takes type (1)
195  // and direction(1)]
196  mesh_pt()->boundary_node_pt(1,j)->pin_position(1,1);
197  mesh_pt()->boundary_node_pt(3,j)->pin_position(1,1);
198 
199  // ...implying that d^2y/ds_0 ds_1=0, too: [pin_position() takes
200  // type (3) and direction (1)]
201  mesh_pt()->boundary_node_pt(1,j)->pin_position(3,1);
202  mesh_pt()->boundary_node_pt(3,j)->pin_position(3,1);
203  }
204  }
205 
206 
207  //Now loop over the sides and apply symmetry conditions
208  //-----------------------------------------------------
209  unsigned n_node_at_side = mesh_pt()->nboundary_node(0);
210  for(unsigned j=0;j<n_node_at_side;j++)
211  {
212 
213  //At the side where theta is 0 (boundary 0), pin in the y displacement
214  //--------------------------------------------------------------------
215 
216  // y=0 for all s_0 (along the tube)...
217  mesh_pt()->boundary_node_pt(0,j)->pin_position(1);
218 
219  // ...implying that dy/ds_0=0 too: [pin_position() takes type (1) and
220  // direction(1)]
221  mesh_pt()->boundary_node_pt(0,j)->pin_position(1,1);
222 
223 
224  //At the side where theta is 0 (boundary 0), apply symm cond for x and z
225  //----------------------------------------------------------------------
226 
227  // d{x,z}/ds_1=0 for all s_0 (along the tube) because of symmetry...
228  // [pin_position() takes type (2) and direction {0,2}]
229  mesh_pt()->boundary_node_pt(0,j)->pin_position(2,0);
230  mesh_pt()->boundary_node_pt(0,j)->pin_position(2,2);
231 
232  // ...implying that d^2{x,z}/ds_0 ds_1=0 too: [pin_position() takes
233  // type (3) and direction {0,2}]
234  mesh_pt()->boundary_node_pt(0,j)->pin_position(3,0);
235  mesh_pt()->boundary_node_pt(0,j)->pin_position(3,2);
236 
237 
238  //At the side where theta is pi/2 (boundary 2), pin in the x displacement
239  //--------------------------------------------------------------------
240 
241  // x=0 for all s_0 (along the tube)...
242  mesh_pt()->boundary_node_pt(2,j)->pin_position(0);
243 
244  // ...implying that dx/ds_0=0 too: [pin_position() takes type (1) and
245  // direction(0)]
246  mesh_pt()->boundary_node_pt(2,j)->pin_position(1,0);
247 
248 
249  //At the side where theta is pi/2 (boundary 2), apply symm cond for y and z
250  //----------------------------------------------------------------------
251 
252  // d{y,z}/ds_1=0 for all s_0 (along the tube) because of symmetry...
253  // [pin_position() takes type (2) and direction {1,2}]
254  mesh_pt()->boundary_node_pt(2,j)->pin_position(2,1);
255  mesh_pt()->boundary_node_pt(2,j)->pin_position(2,2);
256 
257  // ...implying that d^2{y,z}/ds_0 ds_1=0 too: [pin_position() takes
258  // type (3) and direction {1,2}]
259  mesh_pt()->boundary_node_pt(2,j)->pin_position(3,1);
260  mesh_pt()->boundary_node_pt(2,j)->pin_position(3,2);
261 
262  }
263 
264 
265  // Setup displacement control
266  //---------------------------
267 
268  // Choose element in which displacement control is applied: This
269  // one is located about halfway along the tube -- remember that
270  // we've renumbered the elements!
271  unsigned nel_ctrl=0;
272  Vector<double> s_displ_control(2);
273 
274  // Even/odd number of elements in axial direction
275  if (nx%2==1)
276  {
277  nel_ctrl=unsigned(floor(0.5*double(nx))+1.0)*ny-1;
278  s_displ_control[0]=0.0;
279  s_displ_control[1]=1.0;
280  }
281  else
282  {
283  nel_ctrl=unsigned(floor(0.5*double(nx))+1.0)*ny-1;
284  s_displ_control[0]=-1.0;
285  s_displ_control[1]=1.0;
286  }
287 
288  // Controlled element
289  SolidFiniteElement* controlled_element_pt=
290  dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(nel_ctrl));
291 
292 
293  // Fix the displacement in the y (1) direction...
294  unsigned controlled_direction=1;
295 
296  // Pointer to displacement control element
297  DisplacementControlElement* displ_control_el_pt;
298 
299  // Build displacement control element
300  displ_control_el_pt=
301  new DisplacementControlElement(controlled_element_pt,
302  s_displ_control,
303  controlled_direction,
305 
306 
307  // Doc control point
308  Vector<double> xi(2);
309  Vector<double> x(3);
310  controlled_element_pt->interpolated_xi(s_displ_control,xi);
311  controlled_element_pt->interpolated_x(s_displ_control,x);
312  std::cout << std::endl;
313  std::cout << "Controlled element: " << nel_ctrl << std::endl;
314  std::cout << "Displacement control applied at xi = ("
315  << xi[0] << ", " << xi[1] << ")" << std::endl;
316  std::cout << "Corresponding to x = ("
317  << x[0] << ", " << x[1] << ", " << x[2] << ")" << std::endl;
318 
319 
320  // The constructor of the DisplacementControlElement has created
321  // a new Data object whose one-and-only value contains the
322  // adjustable load: Use this Data object in the load function:
323  Global_Physical_Variables::Pext_data_pt=displ_control_el_pt->
324  displacement_control_load_pt();
325 
326  // Add the displacement-control element to the mesh
327  mesh_pt()->add_element_pt(displ_control_el_pt);
328 
329 
330 
331  // Complete build of shell elements
332  //---------------------------------
333 
334  //Find number of shell elements in mesh
335  unsigned n_element = nx*ny;
336 
337  //Explicit pointer to first element in the mesh
338  ELEMENT* first_el_pt = dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(0));
339 
340  //Loop over the elements
341  for(unsigned e=0;e<n_element;e++)
342  {
343  //Cast to a shell element
344  ELEMENT *el_pt = dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(e));
345 
346  //Set the load function
347  el_pt->load_vector_fct_pt() = & Global_Physical_Variables::press_load;
348 
349  //Set the undeformed surface
350  el_pt->undeformed_midplane_pt() = Undeformed_midplane_pt;
351 
352  //The external pressure is external data for all elements
353  el_pt->add_external_data(Global_Physical_Variables::Pext_data_pt);
354 
355  //Pre-compute the second derivatives wrt Lagrangian coordinates
356  //for the first element only
357  if(e==0)
358  {
359  el_pt->pre_compute_d2shape_lagrangian_at_knots();
360  }
361 
362  //Otherwise set the values to be the same as those in the first element
363  //this is OK because the Lagrangian mesh is uniform.
364  else
365  {
366  el_pt->set_dshape_lagrangian_stored_from_element(first_el_pt);
367  }
368  }
369 
370  //Set pointers to two trace nodes, used for output
371  Trace_node_pt = mesh_pt()->finite_element_pt(2*ny-1)->node_pt(3);
372  Trace_node2_pt = mesh_pt()->finite_element_pt(ny)->node_pt(1);
373 
374  // Attach boundary condition elements that enforce clamping condition
375  if (problem_id==0)
376  {
379  }
380 
381 
382  // Do equation numbering
383  cout << std::endl;
384  cout << "# of dofs " << assign_eqn_numbers() << std::endl;
385  cout << std::endl;
386 
387 }
void create_bc_elements(const unsigned &b)
Create clamping bc face elements on the b-th boundary of the Mesh.
Definition: clamped_or_pinned_shell.cc:394
unsigned Nshell
Number of shell elements.
Definition: clamped_or_pinned_shell.cc:125
CircularCylindricalShellMesh< ELEMENT > * mesh_pt()
Overload Access function for the mesh.
Definition: clamped_or_pinned_shell.cc:96
Definition: circular_shell_mesh.template.h:54
virtual double interpolated_x(const Vector< double > &s, const unsigned &i) const
Return FE interpolated coordinate x[i] at local coordinate s.
Definition: elements.cc:3962
virtual double interpolated_xi(const Vector< double > &s, const unsigned &i) const
Definition: elements.cc:7104
EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bfloat16 floor(const bfloat16 &a)
Definition: BFloat16.h:643
list x
Definition: plotDoE.py:28
std::ptrdiff_t j
Definition: tut_arithmetic_redux_minmax.cpp:2

References e(), Eigen::bfloat16_impl::floor(), i, oomph::FiniteElement::interpolated_x(), oomph::SolidFiniteElement::interpolated_xi(), j, Mesh_Parameters::lx, Mesh_Parameters::ly, Mesh_Parameters::nx, Mesh_Parameters::ny, Global_Physical_Variables::Pext_data_pt, Global_Physical_Variables::Prescribed_y, Global_Physical_Variables::press_load(), and plotDoE::x.

◆ ShellProblem() [3/5]

template<class ELEMENT >
ShellProblem< ELEMENT >::ShellProblem ( const unsigned nx,
const unsigned ny,
const double lx,
const double ly 
)

Constructor.

◆ ShellProblem() [4/5]

template<class ELEMENT >
ShellProblem< ELEMENT >::ShellProblem ( const unsigned nx,
const unsigned ny,
const double lx,
const double ly 
)

Constructor.

◆ ~ShellProblem() [2/2]

template<class ELEMENT >
ShellProblem< ELEMENT >::~ShellProblem ( )
inline

Destructor: delete mesh, geometric object.

326  {
327  delete Problem::mesh_pt();
328  delete Undeformed_midplane_pt;
329  }

References ShellProblem< ELEMENT >::Undeformed_midplane_pt.

◆ ShellProblem() [5/5]

template<class ELEMENT >
ShellProblem< ELEMENT >::ShellProblem ( const unsigned nx,
const unsigned ny,
const double lx,
const double ly 
)

Constructor.

Member Function Documentation

◆ actions_after_distribute()

template<class ELEMENT >
void ShellProblem< ELEMENT >::actions_after_distribute ( )
inline

Actions after problem distribution. Need to reset the pointer to stored shape functions on all processors

346  {
347  //If there are any solid elements calculate and store the
348  //shape functions in the first element
349  unsigned n_solid_element = this->solid_mesh_pt()->nelement();
350 
351  //Explicit pointer to first element in the mesh if there is one
352  ELEMENT* first_el_pt =0;
353  if(n_solid_element > 0)
354  {
355  first_el_pt =
356  dynamic_cast<ELEMENT*>(solid_mesh_pt()->element_pt(0));
357  }
358 
359  //Loop over all elements
360  for(unsigned e=0;e<n_solid_element;e++)
361  {
362  ELEMENT* el_pt = dynamic_cast<ELEMENT*>(solid_mesh_pt()->element_pt(e));
363 
364  //Pre-compute the second derivatives wrt Lagrangian coordinates
365  //for the first element only
366  if(e==0)
367  {
368  el_pt->pre_compute_d2shape_lagrangian_at_knots();
369  }
370  //Otherwise set the values to be the same as those in the first element
371  //this is OK because the Lagrangian mesh is uniform.
372  //
373  else
374  {
375  el_pt->set_dshape_lagrangian_stored_from_element(first_el_pt);
376  }
377  }
378  }

References e(), and ShellProblem< ELEMENT >::solid_mesh_pt().

◆ actions_after_newton_solve() [1/5]

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

Actions after solve empty.

Reimplemented from oomph::Problem.

338 {}

◆ actions_after_newton_solve() [2/5]

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

Actions after solve empty.

Reimplemented from oomph::Problem.

102 {}

◆ actions_after_newton_solve() [3/5]

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

Actions after solve empty.

Reimplemented from oomph::Problem.

228 {}

◆ actions_after_newton_solve() [4/5]

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

Actions after solve empty.

Reimplemented from oomph::Problem.

337 {}

◆ actions_after_newton_solve() [5/5]

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

Actions after solve empty.

Reimplemented from oomph::Problem.

107 {}

◆ actions_before_newton_solve() [1/5]

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

Actions before solve empty.

Reimplemented from oomph::Problem.

341 {}

◆ actions_before_newton_solve() [2/5]

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

Actions before solve empty.

Reimplemented from oomph::Problem.

105 {}

◆ actions_before_newton_solve() [3/5]

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

Actions before solve empty.

Reimplemented from oomph::Problem.

231 {}

◆ actions_before_newton_solve() [4/5]

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

Actions before solve empty.

Reimplemented from oomph::Problem.

340 {}

◆ actions_before_newton_solve() [5/5]

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

Actions before solve empty.

Reimplemented from oomph::Problem.

110 {}

◆ create_bc_elements()

template<class ELEMENT >
void ShellProblem< ELEMENT >::create_bc_elements ( const unsigned b)

Create clamping bc face elements on the b-th boundary of the Mesh.

395 {
396  // How many bulk elements are adjacent to boundary b?
397  unsigned n_element = mesh_pt()->nboundary_element(b);
398 
399  std::cout << "Creating " << n_element
400  << " boundary condition elements" << std::endl;
401 
402 
403  // Loop over the bulk elements adjacent to boundary b?
404  for(unsigned e=0;e<n_element;e++)
405  {
406  // Get pointer to the bulk element that is adjacent to boundary b
407  ELEMENT* bulk_elem_pt = dynamic_cast<ELEMENT*>(
408  mesh_pt()->boundary_element_pt(b,e));
409 
410  //What is the index of the face of element e along boundary b
411  int face_index = mesh_pt()->face_index_at_boundary(b,e);
412 
413  // Build the corresponding bc element
414  ClampedHermiteShellBoundaryConditionElement* flux_element_pt = new
415  ClampedHermiteShellBoundaryConditionElement(bulk_elem_pt,face_index);
416 
417  //Add the bc element to the mesh
418  mesh_pt()->add_element_pt(flux_element_pt);
419 
420  } //end of loop over bulk elements adjacent to boundary b
421 
422 }
Scalar * b
Definition: benchVecAdd.cpp:17

References b, and e().

◆ doc_solution()

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

Doc solution.

476 {
477 
478  ofstream some_file;
479  char filename[100];
480 
481  // Loop over shell all elements to get global kinetic and potential energy
482  double global_kin=0;
483  double global_pot=0;
484  double rate_of_work=0;
485  double pot,kin;
486  for (unsigned ielem=0;ielem<Nshell;ielem++)
487  {
488  ELEMENT* el_pt=dynamic_cast<ELEMENT*>(mesh_pt()->element_pt(ielem));
489  el_pt->get_energy(pot,kin);
490  global_kin+=kin;
491  global_pot+=pot;
492  rate_of_work+=el_pt->load_rate_of_work();
493  }
494 
495  //Output the pressure (on the bending scale)
496  trace_file
497  << time_pt()->time() << " "
498  << Trace_node_pt->x(1) << " "
499  << Trace_node2_pt->x(0) << " "
500  << rate_of_work << " "
501  << global_pot << " "
502  << global_kin << " "
503  << global_kin+global_pot << " "
504 // << Global_Physical_Variables::external_pressure()/(pow(0.05,3)/12.0) << " "
505  << Global_Physical_Variables::Pcos/(pow(0.05,3)/12.0) << " "
506  << std::endl;
507 
508  //Output the tube shape
509 
510  sprintf(filename,"%s/shell%i.dat",doc_info.directory().c_str(),
511  doc_info.number());
512  some_file.open(filename);
513  for (unsigned e=0;e<Nshell;e++)
514  {
515  mesh_pt()->finite_element_pt(e)->output(some_file,15);
516  }
517  some_file.close();
518 
519 // //Output the Lagrange multipliers
520 // sprintf(filename,"%s/lagrange_multiplier%i.dat",
521 // doc_info.directory().c_str(),
522 // doc_info.number());
523 // some_file.open(filename);
524 // unsigned n_el=mesh_pt()->nelement();
525 // for (unsigned e=Nshell+1;e<n_el;e++)
526 // {
527 // mesh_pt()->finite_element_pt(e)->output(some_file,15);
528 // }
529 // some_file.close();
530 
531 }
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
double & x(const unsigned &i)
Return the i-th nodal coordinate.
Definition: nodes.h:1060
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
EIGEN_STRONG_INLINE EIGEN_DEVICE_FUNC bfloat16 pow(const bfloat16 &a, const bfloat16 &b)
Definition: BFloat16.h:625
double Pcos
Perturbation pressure.
Definition: steady_ring.cc:56
string filename
Definition: MergeRestartFiles.py:39

References oomph::DocInfo::directory(), e(), MergeRestartFiles::filename, oomph::DocInfo::number(), Global_Physical_Variables::Pcos, and Eigen::bfloat16_impl::pow().

◆ mesh_pt() [1/4]

template<class ELEMENT >
CircularCylindricalShellMesh<ELEMENT>* ShellProblem< ELEMENT >::mesh_pt ( )
inline

Overload Access function for the mesh.

97  {return
98  dynamic_cast<CircularCylindricalShellMesh<ELEMENT>*>(Problem::mesh_pt());
99  }

◆ mesh_pt() [2/4]

template<class ELEMENT >
ShellMesh<ELEMENT>* ShellProblem< ELEMENT >::mesh_pt ( )
inline

Overload Access function for the mesh.

225  {return dynamic_cast<ShellMesh<ELEMENT>*>(Problem::mesh_pt());}

◆ mesh_pt() [3/4]

template<class ELEMENT >
ShellMesh<ELEMENT>* ShellProblem< ELEMENT >::mesh_pt ( )
inline

Overload Access function for the mesh.

334  {return dynamic_cast<ShellMesh<ELEMENT>*>(Problem::mesh_pt());}

◆ mesh_pt() [4/4]

template<class ELEMENT >
CircularCylindricalShellMesh<ELEMENT>* ShellProblem< ELEMENT >::mesh_pt ( )
inline

Overload Access function for the mesh.

102  {return
103  dynamic_cast<CircularCylindricalShellMesh<ELEMENT>*>(Problem::mesh_pt());
104  }

◆ run_it()

template<class ELEMENT >
void ShellProblem< ELEMENT >::run_it

Do parameter study.

Run it...

539 {
540 
541  ofstream some_file;
542  char filename[100];
543 
544  //Increase the maximum number of Newton iterations.
545  //Finding the first buckled solution requires a large(ish) number
546  //of Newton steps -- shells are just a bit twitchy
548  Max_residuals=1.0e6;
549 
550  // Label for output
551  DocInfo doc_info;
552 
553  // Output directory
554  doc_info.set_directory("RESLT");
555 
556  //Open an output trace file
557  sprintf(filename,"%s/trace.dat",doc_info.directory().c_str());
558  ofstream trace_file(filename);
559  trace_file << "VARIABLES=\"time\",\"R_1\",\"R_2\",\"rate of work of load\",\"E_p_o_t\",\"E_k_i_n\",\"E_k_i_n + E_p_o_t\",\"p_c_o_s\"" << std::endl;
560  trace_file << "ZONE" << std::endl;
561 
562 // Set initial timestep
563  double dt=1.0;
564 
565  // Assign impulsive start
567 
568  // Output initial data
569  doc_solution(doc_info,trace_file);
570 
571 
572 
573 
574  // Reduce number of timesteps for validation
575  unsigned nstep=200;
576  if (CommandLineArgs::Argc>1) nstep=4;
577 
578 
579  // Time integration loop
580  for(unsigned i=1;i<=nstep;i++)
581  {
582 
583  // Switch off perturbation pressure
585  {
586  // Perturbation pressure
588  }
589 
590  // Solve
592 
593  // Doc solution
594  doc_info.number()++;
595  doc_solution(doc_info,trace_file);
596 
597  }
598 
599  //Close the trace file
600  trace_file.close();
601 
602 }
void doc_solution(DocInfo &doc_info, ofstream &trace_file)
Doc solution.
Definition: oscillating_shell.cc:474
Definition: oomph_utilities.h:499
void set_directory(const std::string &directory)
Definition: oomph_utilities.cc:298
unsigned Max_newton_iterations
Maximum number of Newton iterations.
Definition: problem.h:599
void assign_initial_values_impulsive()
Definition: problem.cc:11499
double Max_residuals
Definition: problem.h:610
void unsteady_newton_solve(const double &dt)
Definition: problem.cc:10953
double & time()
Return the current value of continuous time.
Definition: problem.cc:11531
double T_pcos_end
Only keep p_cos switched on until T_pcos_end.
Definition: oscillating_shell.cc:58
int Argc
Number of arguments + 1.
Definition: oomph_utilities.cc:407

References oomph::CommandLineArgs::Argc, oomph::DocInfo::directory(), MergeRestartFiles::filename, i, oomph::Locate_zeta_helpers::Max_newton_iterations, oomph::DocInfo::number(), Global_Physical_Variables::Pcos, oomph::DocInfo::set_directory(), and Global_Physical_Variables::T_pcos_end.

◆ solid_mesh_pt()

template<class ELEMENT >
ShellMesh<ELEMENT>* ShellProblem< ELEMENT >::solid_mesh_pt ( )
inline

Overload Access function for the mesh.

335  {return dynamic_cast<ShellMesh<ELEMENT>*>(this->Solid_mesh_pt);}

References ShellProblem< ELEMENT >::Solid_mesh_pt.

Referenced by ShellProblem< ELEMENT >::actions_after_distribute().

◆ solve() [1/4]

template<class ELEMENT >
void ShellProblem< ELEMENT >::solve
628 {
629 
630  //Increase the maximum number of Newton iterations.
631  //Finding the first buckled solution requires a large(ish) number
632  //of Newton steps -- shells are just a bit twitchy
634 
635  //Open an output trace file
636  char filename[100];
637  sprintf(filename,"trace_on_proc%i.dat",this->communicator_pt()->my_rank());
638  ofstream trace(filename);
639 
640  //Change in displacemenet
641  double disp_incr = 0.05;
642 
643  //Gradually compress the tube by decreasing the value of the prescribed
644  //position from 1 to zero in steps of 0.05 initially and then 0.1
645  for(unsigned i=1;i<13;i++)
646  {
647  //By the time we reach the second time round increase the incremenet
648  if(i==3) {disp_incr = 0.1;}
649  //Reduce prescribed y by our chosen increment
651 
652  cout << std::endl << "Increasing displacement: Prescribed_y is "
654 
655  // Solve
656  newton_solve();
657 
658  //Output the pressure
659  trace << Global_Physical_Variables::external_pressure()/(pow(0.05,3)/12.0)
660  << " "
661  //Position of first trace node
662  << Trace_node_pt->x(0) << " " << Trace_node_pt->x(1) << " "
663  //Position of second trace node
664  << Trace_node2_pt->x(0) << " " << Trace_node2_pt->x(1) << std::endl;
665  }
666 
667  //Close the trace file
668  trace.close();
669 
670  //Output the tube shape in the most strongly collapsed configuration
671  sprintf(filename,"final_shape_on_proc%i.dat",
672  this->communicator_pt()->my_rank());
673  ofstream file(filename);
674  solid_mesh_pt()->output(file,5);
675  file.close();
676 
677 
678  //Switch from displacement control to arc-length continuation and
679  //trace back up the solution branch
680 
681  //Now pin the external pressure
683 
684  //Re-assign the equation numbers
685  cout << std::endl;
686  cout << "-----------------ARC-LENGTH CONTINUATION --------------"
687  << std::endl;
688  cout << "# of dofs " << assign_eqn_numbers() << std::endl;
689  cout << "-------------------------------------------------------"
690  << std::endl;
691  cout << std::endl;
692 
693  //Set the maximum number of Newton iterations to something more reasonable
695 
696  //Set the desired number of Newton iterations per arc-length step
698 
699  //Set the proportion of the arc length
701 
702  //Set an initial value for the step size
703  double ds = -0.5;
704 
705  //Open a different trace file
706  sprintf(filename,"trace_disp_on_proc%i.dat",
707  this->communicator_pt()->my_rank());
708  trace.open(filename);
709  //Take fifteen continuation steps
710  for(unsigned i=0;i<15;i++)
711  {
714 
715  //Output the pressure
716  trace << Global_Physical_Variables::external_pressure()/(pow(0.05,3)/12.0)
717  << " "
718  //Position of first trace node
719  << Trace_node_pt->x(0) << " " << Trace_node_pt->x(1) << " "
720  //Position of second trace node
721  << Trace_node2_pt->x(0) << " " << Trace_node2_pt->x(1) << std::endl;
722  }
723 
724  //Close the trace file
725  trace.close();
726 
727 }
void pin(const unsigned &i)
Pin the i-th stored variable.
Definition: nodes.h:385
double Desired_proportion_of_arc_length
Proportion of the arc-length to taken by the parameter.
Definition: problem.h:734
OomphCommunicator * communicator_pt()
access function to the oomph-lib communicator
Definition: problem.h:1246
void newton_solve()
Use Newton method to solve the problem.
Definition: problem.cc:8783
unsigned Desired_newton_iterations_ds
Definition: problem.h:784
double arc_length_step_solve(double *const &parameter_pt, const double &ds, const unsigned &max_adapt=0)
Definition: problem.cc:10294
double & external_pressure()
Access function to value of external pressure.
Definition: steady_ring.cc:79

References Global_Physical_Variables::external_pressure(), MergeRestartFiles::filename, i, oomph::Locate_zeta_helpers::Max_newton_iterations, Global_Physical_Variables::Pext_data_pt, oomph::Data::pin(), Eigen::ArrayBase< Derived >::pow(), and Global_Physical_Variables::Prescribed_y.

◆ solve() [2/4]

template<class ELEMENT >
void ShellProblem< ELEMENT >::solve ( )

◆ solve() [3/4]

template<class ELEMENT >
void ShellProblem< ELEMENT >::solve ( )

◆ solve() [4/4]

template<class ELEMENT >
void ShellProblem< ELEMENT >::solve ( string &  dir_name)

Do parameter study.

430 {
431  ofstream some_file;
432  char filename[100];
433 
434  //Increase the maximum number of Newton iterations.
435  //Finding the first buckled solution requires a large(ish) number
436  //of Newton steps -- shells are just a bit twitchy
438  Max_residuals=1.0e6;
439 
440  // Label for output
441  DocInfo doc_info;
442 
443  // Output directory
444  doc_info.set_directory(dir_name);
445 
446  //Open an output trace file
447  sprintf(filename,"%s/trace.dat",doc_info.directory().c_str());
448  ofstream trace_file(filename);
449  trace_file << "VARIABLES=\"p_e_x_t\",\"R_1\",\"R_2\"" << std::endl;
450  trace_file << "ZONE" << std::endl;
451 
452  //Gradually compress the tube by decreasing the value of the prescribed
453  //position
454  for(unsigned i=1;i<11;i++)
455  {
456 
458 
459  // Solve
460  newton_solve();
461 
462  //Output the pressure (on the bending scale)
463  trace_file
464  << Global_Physical_Variables::external_pressure()/(pow(0.05,3)/12.0) << " "
465  << Trace_node_pt->x(1) << " "
466  << Trace_node2_pt->x(0) << std::endl;
467 
468 
469  //Output the tube shape
470  sprintf(filename,"%s/shell%i.dat",doc_info.directory().c_str(),
471  doc_info.number());
472  some_file.open(filename);
473  for (unsigned e=0;e<Nshell;e++)
474  {
475  mesh_pt()->finite_element_pt(e)->output(some_file,15);
476  }
477  some_file.close();
478 
479  //Output the Lagrange multipliers
480  sprintf(filename,"%s/lagrange_multiplier%i.dat",
481  doc_info.directory().c_str(),
482  doc_info.number());
483  some_file.open(filename);
484  unsigned n_el=mesh_pt()->nelement();
485  for (unsigned e=Nshell+1;e<n_el;e++)
486  {
487  mesh_pt()->finite_element_pt(e)->output(some_file,15);
488  }
489  some_file.close();
490 
491  // Increment counter for output files
492  doc_info.number()++;
493 
494  // Reset perturbation
496  }
497 
498  //Close the trace file
499  trace_file.close();
500 
501 }

References oomph::DocInfo::directory(), e(), Global_Physical_Variables::external_pressure(), MergeRestartFiles::filename, i, oomph::Locate_zeta_helpers::Max_newton_iterations, oomph::DocInfo::number(), Global_Physical_Variables::Pcos, Eigen::bfloat16_impl::pow(), Global_Physical_Variables::Prescribed_y, and oomph::DocInfo::set_directory().

Member Data Documentation

◆ Displacement_control_mesh_pt

template<class ELEMENT >
Mesh* ShellProblem< ELEMENT >::Displacement_control_mesh_pt
private

Pointer to the mesh that contains the displacement control element.

Referenced by ShellProblem< ELEMENT >::~ShellProblem().

◆ Nshell

template<class ELEMENT >
unsigned ShellProblem< ELEMENT >::Nshell
private

Number of shell elements.

◆ Solid_mesh_pt

template<class ELEMENT >
Mesh* ShellProblem< ELEMENT >::Solid_mesh_pt
private

◆ Trace_node2_pt

template<class ELEMENT >
Node * ShellProblem< ELEMENT >::Trace_node2_pt
private

Second trace node.

◆ Trace_node_pt

template<class ELEMENT >
Node * ShellProblem< ELEMENT >::Trace_node_pt
private

First trace node.

◆ Undeformed_midplane_pt

template<class ELEMENT >
GeomObject * ShellProblem< ELEMENT >::Undeformed_midplane_pt
private

Pointer to GeomObject that specifies the undeformed midplane.

Referenced by ShellProblem< ELEMENT >::~ShellProblem().


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