interaction/pseudo_solid_collapsible_tube/pseudo_solid_collapsible_tube.cc File Reference

#include "generic.h"
#include "solid.h"
#include "constitutive.h"
#include "navier_stokes.h"
#include "multi_physics.h"
#include "meshes/simple_cubic_mesh.h"
#include "meshes/quarter_tube_mesh.h"

Classes
class	oomph::PseudoElasticBulkElement< ELEMENT >
class	oomph::FaceGeometry< PseudoElasticBulkElement< ELEMENT > >
	FaceGeometry of wrapped element is the same as the underlying element. More...
class	WallMesh< ELEMENT >
class	FluidMesh< ELEMENT >
class	PseudoElasticCollapsibleChannelProblem< FLUID_ELEMENT, SOLID_ELEMENT >

Namespaces
	oomph
	DRAIG: Change all instances of (SPATIAL_DIM) to (DIM-1).
	Global_Parameters
	Namespace for global parameters.

Typedefs
typedef PseudoElasticCollapsibleChannelProblem< RefineablePseudoSolidNodeUpdateElement< RefineableQTaylorHoodElement< 3 >, PseudoElasticBulkElement< RefineableQPVDElement< 3, 3 > > >, RefineableQPVDElement< 3, 3 > >	PseudoElasticFSIProblem

Functions
void	Global_Parameters::press_load (const Vector< double > &xi, const Vector< double > &x, const Vector< double > &N, Vector< double > &load)
	Load function for wall. More...
void	set_pseudo_elastic_fsi_solver (PseudoElasticFSIProblem &problem)
int	main (int argc, char *argv[])

Variables
double	Global_Parameters::L =12.0
	tube length More...
double	Global_Parameters::Lup =1.5
	upstream length More...
double	Global_Parameters::Ldown =3.0
	downstream length More...
double	Global_Parameters::H =0.15
	wall thickness More...
unsigned	Global_Parameters::N_slice =4
	number of axial slices in fluid mesh More...
ConstitutiveLaw *	Global_Parameters::Constitutive_law_wall_pt =0
	Pointer to constitutive law for the wall. More...
ConstitutiveLaw *	Global_Parameters::Constitutive_law_pseudo_elastic_pt =0
	Pointer to constitutive law for the pseudo elastic node update elements. More...
double	Global_Parameters::Nu_wall =0.3
double	Global_Parameters::Nu_pseudo_elastic =0.1
double	Global_Parameters::Pcos = 2.0e-4
	Perturbation pressure. More...
double	Global_Parameters::Lambda_sq =0.0
	Timescale ratio (non-dim density) for solid. More...
double	Global_Parameters::T =1.0
	Period of periodic variation in inflow pressure. More...
double	Global_Parameters::Dt =0.01
	Timestep. More...
unsigned	Global_Parameters::Nstep_per_period =40
	Number of steps per period. More...
unsigned	Global_Parameters::Nperiod =5
	Number of periods. More...

Typedef Documentation

PseudoElasticFSIProblem

typedef PseudoElasticCollapsibleChannelProblem< RefineablePseudoSolidNodeUpdateElement<RefineableQTaylorHoodElement<3>, PseudoElasticBulkElement<RefineableQPVDElement<3,3> > >,RefineableQPVDElement<3,3> > PseudoElasticFSIProblem

Function Documentation

main()

int main	(	int argc	,
		char *	argv[]
	)

Demonstrate how to solve a 3D FSI problem with pseudo-solid node update for the fluid mesh.

{       
 // Store command line arguments
 CommandLineArgs::setup(argc,argv);
 
 // Create generalised Hookean constitutive equations
 Global_Parameters::Constitutive_law_wall_pt = 
  new GeneralisedHookean(&Global_Parameters::Nu_wall);
 
Global_Parameters::Constitutive_law_pseudo_elastic_pt = 
  new GeneralisedHookean(&Global_Parameters::Nu_pseudo_elastic);
 
 //Set up the problem
 PseudoElasticFSIProblem problem; 
 
 // Doc info
 DocInfo doc_info;
 doc_info.set_directory("RESLT");
 
 // Doc initial configuration
 problem.doc_solution(doc_info);
 oomph_info << "Before manual refine: doc_info.number()=" 
            << doc_info.number() << std::endl;
 doc_info.number()++;
 
 // Do manual non-uniform refinement
 if (CommandLineArgs::Argc!=1)
  {
   problem.adapt(); 
  }
 
 // Doc afterwards
 problem.doc_solution(doc_info);
 oomph_info << "After manual refine: doc_info.number()=" 
            << doc_info.number() << std::endl;
 doc_info.number()++;
 
 // Setup timestepping
 double dt = Global_Parameters::T/double(Global_Parameters::Nstep_per_period);
 Global_Parameters::Dt=dt;
 problem.initialise_dt(dt);
 problem.set_initial_condition();
 
  // set
 if (CommandLineArgs::Argc>2) 
  {
   set_pseudo_elastic_fsi_solver(problem);
  }
 
 // Steady run
 problem.steady_run(doc_info);
 
 // Unteady run
 problem.unsteady_run(doc_info);
 
 // clean up
 if (CommandLineArgs::Argc>2) 
  {
   IterativeLinearSolver* solver_pt = dynamic_cast<IterativeLinearSolver*>(problem.linear_solver_pt());
   delete solver_pt->preconditioner_pt();
   delete solver_pt;
  }
} // end_of_main

References oomph::CommandLineArgs::Argc, Global_Parameters::Constitutive_law_pseudo_elastic_pt, Global_Parameters::Constitutive_law_wall_pt, Global_Parameters::Dt, Global_Parameters::Nstep_per_period, Global_Parameters::Nu_pseudo_elastic, Global_Parameters::Nu_wall, oomph::DocInfo::number(), oomph::oomph_info, oomph::IterativeLinearSolver::preconditioner_pt(), problem, oomph::DocInfo::set_directory(), set_pseudo_elastic_fsi_solver(), Flag_definition::setup(), and Global_Parameters::T.

set_pseudo_elastic_fsi_solver()

void set_pseudo_elastic_fsi_solver

(

PseudoElasticFSIProblem &

problem

)

helper method to return GMRES preconditioned with the pseudo-elastic FSI preconditioner

{
// setup the solver
//#ifdef OOMPH_HAS_TRILINOS
// TrilinosAztecOOSolver* solver_pt = new
//  TrilinosAztecOOSolver;
// solver_pt->solver_type() = TrilinosAztecOOSolver::GMRES;
//#else
 GMRES<CRDoubleMatrix>* solver_pt = new GMRES<CRDoubleMatrix>;
//#endif
 solver_pt->tolerance()=1e-8;
 
 // preconditioner
 PseudoElasticFSIPreconditioner* prec_pt = new
  PseudoElasticFSIPreconditioner(3, &problem);
 
 // meshes
 prec_pt->set_fluid_and_pseudo_elastic_mesh_pt(problem.fluid_mesh_pt());
 prec_pt->set_solid_mesh_pt(problem.wall_mesh_pt());
 prec_pt->set_lagrange_multiplier_mesh_pt(
  problem.lagrange_multiplier_mesh_pt());
 
/*
 // inexact pseudo-solid preconditioning
 prec_pt->pseudo_elastic_preconditioner_pt()->elastic_preconditioner_type()
  = PseudoElasticPreconditioner::
  Block_upper_triangular_preconditioner;
 prec_pt->pseudo_elastic_preconditioner_pt()
  ->set_elastic_subsidiary_preconditioner
  (Pseudo_Elastic_Preconditioner_Subsidiary_Operator_Helper
   ::get_elastic_preconditioner);
 prec_pt->pseudo_elastic_preconditioner_pt()
  ->set_lagrange_multiplier_subsidiary_preconditioner
  (Pseudo_Elastic_Preconditioner_Subsidiary_Operator_Helper
   ::get_lagrange_multiplier_preconditioner);
 
 // inexact "real" solid preconditioning
 BlockTriangularPreconditioner<CRDoubleMatrix>*
  solid_prec_pt = new BlockTriangularPreconditioner<CRDoubleMatrix>;
 solid_prec_pt->set_subsidiary_preconditioner_function
  (Real_Solid_Preconditioner_Helper::get_preconditioner);
 
 // inexact navier stokes preconditioning
  NavierStokesSchurComplementPreconditioner*
  ns_prec_pt = prec_pt->navier_stokes_schur_complement_preconditioner_pt();
 prec_pt->enable_navier_stokes_schur_complement_preconditioner();
 
 // ns momentum
 BlockDiagonalPreconditioner<CRDoubleMatrix>*
  f_prec_pt = new BlockDiagonalPreconditioner<CRDoubleMatrix>;
 f_prec_pt->set_subsidiary_preconditioner_function
  (LSC_Preconditioner_Helper::set_hypre_preconditioner);
 ns_prec_pt->set_f_preconditioner(f_prec_pt);
 
 // ns pressure poisson
 HyprePreconditioner* p_prec_pt = new HyprePreconditioner;
 p_prec_pt->set_amg_iterations(2);
 p_prec_pt->amg_using_simple_smoothing();
 p_prec_pt->amg_simple_smoother() = 3;
 p_prec_pt->hypre_method() = HyprePreconditioner::BoomerAMG;
 p_prec_pt->amg_strength() = 0.25;
 // Use 6 for parallel (Falgout) or 3 (Ruge Stuben)or serial
 p_prec_pt->amg_coarsening() = 6; 
 p_prec_pt->disable_doc_time();
 ns_prec_pt->set_p_preconditioner(p_prec_pt);
*/
 
 // and pass prec to solver and solver to problem
 solver_pt->preconditioner_pt() = prec_pt;
 problem.linear_solver_pt() = solver_pt;
}

References e(), oomph::IterativeLinearSolver::preconditioner_pt(), problem, oomph::PseudoElasticFSIPreconditioner::set_fluid_and_pseudo_elastic_mesh_pt(), oomph::PseudoElasticFSIPreconditioner::set_lagrange_multiplier_mesh_pt(), oomph::PseudoElasticFSIPreconditioner::set_solid_mesh_pt(), and oomph::IterativeLinearSolver::tolerance().

Referenced by main().