====== poincare ====== Compute crossings of a Poincare section in a pre-computed trajectory. The Poincare condition is %%(u(t) - ueqb, e) == 0%%, i.e. the displacement from a given point %%ueqb%% is orthogonal to a direction %%e%%. The program checks for changes in sign in the inner product in the stored trajectory (with large %%dT%%, and when it finds one, goes back and reintegrates with small %%dt%%, and then does quadratic interpolation when the fine-scale trajectory recrosses the section. ===== usage ===== poincare -T0 0 -T1 400 -d datadir ueqb.ff e.ff reads a trajectory u0, u1, ..., u400 stored in %%datadir%%, computes the crossings of the Poincare section %%(u(t)-ueqb, e) == 0%%, and saves them in a %%poincare/%% directory. ===== options ===== options : -T0 --T0 default == 0 start time -T1 --T1 default == 100 end time -dT --dT default == 1 save interval -d --datadir default == data/ flowfield series directory -o --poindir default == poincare/ output dir for Poincare crossings -vdt --variabledt adjust dt to keep CFLmin<=CFL default == 0.03125 timestep -dtmin --dtmin default == 0.001 minimum time step -dtmax --dtmax default == 0.05 maximum time step -CFLmin --CFLmin default == 0.4 minimum CFL number -CFLmax --CFLmax default == 0.6 maximum CFL number -ts --timestepping default == sbdf3 timestepping algorithm -nl --nonlinearity default == rot method of calculating nonlinearity -R --Reynolds default == 400 Reynolds number -c --channel channelflow instead of plane Couette -b --bulkvelocity hold bulk velocity fixed, not pressure gradient -P --dPdx default == 0 value for fixed pressure gradient -U --Ubulk default == 0 value for fixed bulk velocity (trailing arg 2) fixed pt of poincare map (trailing arg 1) e, for poincare condition (u(t)-u*, e) = 0 ===== Usage example ===== Here is the basic idea of usage. Suppose you have an equilibrium ueqb with an unstable complex eigenvalue pair whose eigenfunctions are ef1 and ef2. Integrate a trajectory with initial condition %%u(0) = ueqb + 0.001 ef1%% as follows addfields 1 ueqb 0.001 ef1 ueqb_001ef1 couette -T0 0 -T1 1000 -o data-ueqb_0001ef1 ueqb_001ef1 Then you can define a Poincare section as, say, %%(u-ueqb, ef2) = 0%%, that is, the displacement from the equilibrium is orthogonal to the ef2 eigenfunction. To compute crossings of the Poincare section within the previously integrated trajectory, run poincare -T0 0 -T1 1000 -d data-ueqb_0001ef1 ueqb ef2 ===== Hiccups ===== :-) Hey John, I was wondering if there was any particular reason you put the cfpause() line in your code other than it was helpful in writing the code. --- //[[dustin.spieker@gatech.edu|Dustin Spieker]] 2009-03-02 10:24//