adaptive_hopf.cc File Reference
#include <iostream>
#include <fstream>
#include <cstdio>
#include "generic.h"
#include "navier_stokes.h"
#include "assert.h"

Classes

class  GeneralEllipse
 My own Ellipse class. More...
 
class  RectangleWithHoleDomain
 Rectangular domain with circular whole. More...
 
class  RectangleWithHoleMesh< ELEMENT >
 Domain-based mesh for rectangular mesh with circular hole. More...
 
class  RefineableRectangleWithHoleMesh< ELEMENT >
 
class  FlowAroundCylinderProblem< ELEMENT >
 Flow around a cylinder in rectangular domain. More...
 

Namespaces

 Global_Parameters
 Namespace for global parameters.
 

Functions

int main ()
 Driver. More...
 

Variables

double Global_Parameters::Re =75.0
 reynolds number More...
 
double Global_Parameters::B =0.7
 
double Global_Parameters::Alpha =0.0
 Peakiness parameter for pressure load. More...
 
bool Global_Parameters::Read_in_eigenfunction_from_disk = true
 

Function Documentation

◆ main()

int main ( )

Driver.

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

1248 {
1249  // Length and height of domain
1250  double length = 14.0;
1251  double height = 1.0/Global_Parameters::B;
1252 
1253  //Create a new ellipse object with equal semi-minor and semi-major axes
1254  //to be the central cylinder
1255  GeneralEllipse* cylinder_pt =
1256  new GeneralEllipse(4.5,0.0,0.5,0.5);
1257 
1258  // Create Problem
1260  <RefineableQCrouzeixRaviartElement<2> > problem(cylinder_pt,length,height);
1261 
1262  //Refine the problem uniformly a couple of times
1263  problem.refine_uniformly(); problem.refine_uniformly();
1264 
1265  // Solve adaptively with up to two rounds of refinement
1266  problem.newton_solve(2);
1267 
1268  //Now increase the Reynolds number to 90 in steps of 5
1269  //without any further adaptation
1270  for(unsigned i=0;i<3;i++)
1271  {
1272  Global_Parameters::Re += 5.0;
1273  problem.newton_solve();
1274  }
1275 
1276  //Assign memory for the eigenvalues and eigenvectors
1277  Vector<DoubleVector> eigenvectors;
1278  double frequency = 0.0;
1279 
1280  //If we are reading in from the disk
1282  {
1283  //Read in the eigenvalue from the data file
1284  std::ifstream input("eigen.dat");
1285  input >> frequency;
1286 
1287  //Read in the eigenvector from the data file
1288  const unsigned n_dof = problem.ndof();
1289  eigenvectors.resize(2);
1290  LinearAlgebraDistribution dist(problem.communicator_pt(),n_dof,false);
1291  //Rebuild the vector
1292  eigenvectors[0].build(&dist,0.0);
1293  eigenvectors[1].build(&dist,0.0);
1294 
1295  for(unsigned n=0;n<n_dof;n++)
1296  {
1297  input >> eigenvectors[0][n];
1298  input >> eigenvectors[1][n];
1299  }
1300  input.close();
1301  }
1302  //Otherwise solve the eigenproblem
1303  else
1304  {
1305  Vector<std::complex<double> > eigenvalues;
1306  //Now solve the eigenproblem
1307  problem.solve_eigenproblem(6,eigenvalues,eigenvectors);
1308  frequency = eigenvalues[0].imag();
1309  }
1310 
1311  //Try to find the Hopf exactly by using the initial
1312  //guesses for the eigenvalue and eigenvector from
1313  //the data file
1314  problem.activate_hopf_tracking(&Global_Parameters::Re,
1315  frequency,
1316  eigenvectors[0],
1317  eigenvectors[1]);
1318  //Solve the problem
1319  problem.newton_solve();
1320  //Report the value of the bifurcation
1321  std::cout << "Hopf bifurcation found at " << Global_Parameters::Re << " " <<
1322  problem.dof(problem.ndof()-1) << std::endl;
1323 
1324  //Solve it again with one round of adaptation
1325  problem.newton_solve(1);
1326 
1327  std::cout << "Hopf bifurcation found at " << Global_Parameters::Re << " " <<
1328  problem.dof(problem.ndof()-1) << std::endl;
1329  problem.mesh_pt()->output("bif_soln.dat");
1330 
1331  //Output the value of the critical Reynolds number into a data file
1332  std::ofstream trace("trace.dat");
1333  trace << Global_Parameters::Re << " " <<
1334  problem.dof(problem.ndof()-1) << "\n";
1335  trace.close();
1336 }
int i
Definition: BiCGSTAB_step_by_step.cpp:9
const unsigned n
Definition: CG3DPackingUnitTest.cpp:11
Flow around a cylinder in rectangular domain.
Definition: adaptive_hopf.cc:1042
My own Ellipse class.
Definition: adaptive_hopf.cc:80
Definition: linear_algebra_distribution.h:64
Definition: oomph-lib/src/generic/Vector.h:58
bool Read_in_eigenfunction_from_disk
Definition: adaptive_hopf.cc:65
double Re
reynolds number
Definition: adaptive_hopf.cc:54
double B
Definition: adaptive_hopf.cc:57
double height(const double &x)
Height of domain.
Definition: simple_spine_channel.cc:429
Constructor for SteadyAxisymAdvectionDiffusion problem
Definition: steady_axisym_advection_diffusion.cc:213

References Global_Parameters::B, Global_Physical_Variables::height(), i, n, problem, Global_Parameters::Re, and Global_Parameters::Read_in_eigenfunction_from_disk.