fsi_pseudo_solid_collapsible_channel_adapt.cc File Reference

#include "generic.h"
#include "navier_stokes.h"
#include "beam.h"
#include "solid.h"
#include "constitutive.h"
#include "meshes/one_d_lagrangian_mesh.h"
#include "meshes/collapsible_channel_mesh.h"

Classes
class	ElasticRefineableCollapsibleChannelMesh< ELEMENT >
	Upgrade mesh to solid mesh. More...
class	UndeformedWall
class	FSICollapsibleChannelProblem< ELEMENT >
	Problem class. More...

Namespaces
	BL_Squash
	Global_Physical_Variables
	Global variables.

Functions
double	BL_Squash::squash_fct (const double &s)
void	Global_Physical_Variables::prescribed_traction (const double &t, const Vector< double > &x, const Vector< double > &n, Vector< double > &traction)
	Traction applied on the fluid at the left (inflow) boundary. More...
void	Global_Physical_Variables::load (const Vector< double > &xi, const Vector< double > &x, const Vector< double > &N, Vector< double > &load)
	Load function: Apply a constant external pressure to the beam. More...
int	main (int argc, char *argv[])

Function Documentation

main()

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

Driver code for a collapsible channel problem with FSI. Presence of command line arguments indicates validation run with coarse resolution and small number of timesteps.

{
 // Store command line arguments
 CommandLineArgs::setup(argc,argv);
  
 // Reduction in resolution for validation run?
 unsigned coarsening_factor=4; 
 if (CommandLineArgs::Argc>1)
  {
   coarsening_factor=4;
  }
 
 // Number of elements in the domain
 unsigned nup=20/coarsening_factor;
 unsigned ncollapsible=40/coarsening_factor;
 unsigned ndown=40/coarsening_factor;
 unsigned ny=16/coarsening_factor;
 
 // Length of the domain
 double lup=5.0;
 double lcollapsible=10.0;
 double ldown=10.0;
 double ly=1.0;
 
 // Set external pressure (on the wall stiffness scale). 
 Global_Physical_Variables::P_ext = 1.0e-1;
 
 // Pressure on the left boundary: This is consistent with steady
 // Poiseuille flow
 Global_Physical_Variables::P_up=12.0*(lup+lcollapsible+ldown);
 
 
#ifdef TAYLOR_HOOD
 
 // Build the problem with QtaylorHoodElements
 FSICollapsibleChannelProblem
  <RefineablePseudoSolidNodeUpdateElement<
  RefineableQTaylorHoodElement<2>,
  RefineableQPVDElement<2,3> > >
  problem(nup, ncollapsible, ndown, ny, 
          lup, lcollapsible, ldown, ly);
 
#else
 
 // Build the problem with QCrouzeixRaviartElements
 FSICollapsibleChannelProblem
  <RefineablePseudoSolidNodeUpdateElement<
  RefineableQCrouzeixRaviartElement<2>,
  RefineableQPVDElement<2,3> > >
  problem(nup, ncollapsible, ndown, ny, 
          lup, lcollapsible, ldown, ly);
 
#endif
 
 // Timestep. Note: Preliminary runs indicate that the period of
 // the oscillation is about 1 so this gives us 40 steps per period.
 double dt=1.0/40.0; 
 
 // Initial time for the simulation
 double t_min=0.0;
 
 // Maximum time for simulation
 double t_max=3.5; 
 
 // Initialise timestep 
 problem.time_pt()->time()=t_min;
 problem.initialise_dt(dt);
 
 // Apply initial condition
 problem.set_initial_condition();
 
 //Set output directory
 DocInfo doc_info;
 doc_info.set_directory("RESLT");
 
 // Open a trace file 
 ofstream trace_file;
 char filename[100];   
 sprintf(filename,"%s/trace.dat",doc_info.directory().c_str());
 trace_file.open(filename);
 
 // Output the initial condition
 problem.doc_solution(doc_info, trace_file);
 
 // Increment step number
 doc_info.number()++;
 
 // Find number of timesteps (reduced for validation)
 unsigned nstep = unsigned((t_max-t_min)/dt);
 if (CommandLineArgs::Argc>1)
  {
   nstep=3;
  }
 
 
 // Set targets for spatial adaptivity
 problem.bulk_mesh_pt()->max_permitted_error()=1.0e-3;
 problem.bulk_mesh_pt()->min_permitted_error()=1.0e-5;
 
 // Overwrite for validation run
 if (CommandLineArgs::Argc>1)
  {
   // Set targets for spatial adaptivity
   problem.bulk_mesh_pt()->max_permitted_error()=0.5e-2;
   problem.bulk_mesh_pt()->min_permitted_error()=0.5e-4;
  }
 
 // When performing the first timestep, we can adapt the mesh as many times
 // as we want because the initial condition can be re-set
 unsigned max_adapt=3; 
 bool first=true;
 
 // Timestepping loop
 for (unsigned istep=0;istep<nstep;istep++)
  {
   // Solve the problem
   problem.unsteady_newton_solve(dt,max_adapt,first);
   
   // Outpt the solution
   problem.doc_solution(doc_info, trace_file);
   
   // Step number
   doc_info.number()++;
 
   // We've done the first step
   first=false;
   max_adapt=1; 
  }
 
 // Close trace file.
 trace_file.close();
 
}//end of main