findbif

 find bifpoint starting at p such that ineg=p.ineg+ichange.
 ichange=1 or -1
 note: no "fold detection" since ultimate bifdetec via cont
       still only works for simple Evals crossing imag axis
 ideally, user should already which direction the index change is

0001 function p=findbif(p,ichange)
0002 % find bifpoint starting at p such that ineg=p.ineg+ichange.
0003 % ichange=1 or -1
0004 % note: no "fold detection" since ultimate bifdetec via cont
0005 %       still only works for simple Evals crossing imag axis
0006 % ideally, user should already which direction the index change is
0007 
0008 smod=p.smod;      % saving parameters, which we will change in this
0009 bif=p.bifchecksw; % function to set them back at the end;
0010 ds=p.ds; pnamesw=p.pnamesw; headfu=p.headfu; ufu=p.ufu; count=p.count; 
0011 timesw=p.timesw; pmod=p.pmod; nsteps=p.nsteps; dsmin=p.dsmin; 
0012 
0013 p.smod=0; % here we do not want to save any point,
0014 % because we are only looking for the bifpoint
0015 p.bifchecksw=0; % the bifpoint will we calc. at the end
0016 p.pnamesw=0;    % no savings
0017 p.headfu=@findbifheadfu; % no outputs
0018 p.ufu=@findbiffu; p.timesw=0; p.pmod=0; p.nsteps=1; 
0019 p.dsmin=p.dsmin/10; % might need very small stepsize!
0020 
0021 st=0; % while loop stop parameter
0022 af=0; % already found a solution with |index-change|>1
0023 ineg1=p.ineg; % number of negative eigenvalues of starting point
0024 while (st==0 && abs(p.ds)>p.dsmin) 
0025   ft=p; ft=cont(ft);
0026   figure(ft.ifig); plot(real(ft.muv),imag(ft.muv),'*'); 
0027   if(p.vsw>0)    
0028      disp(['lam=' num2str(ft.lam) ' # negative EVals=' num2str(ft.ineg)] ); 
0029   end
0030   if ft.ineg==ineg1 % index didn't change, next step!
0031      p=ft; 
0032      if af==1; p.ds=p.ds/4;end
0033      af=0; % assume that next solution will have no large index-change
0034   end
0035   if abs(ft.ineg-ineg1)>1 % gone too far
0036      p.ds=p.ds/2; af=1; end
0037   if ft.ineg==ineg1+ichange; % point found!
0038      st=1; end     
0039   if (ft.ineg-ineg1)*ichange<0; % gone into wrong direction, flip!
0040       p.ds=-p.ds; 
0041   end
0042 end   % while
0043 % now bifurcation (unless index change due to fold) should be found
0044 % with current p.ds, first restore parameters
0045 p.smod=smod; p.bifchecksw=1; p.headfu=headfu; 
0046 p.ufu=ufu; p.timesw=timesw; p.pmod=pmod; p.count=count;
0047 p=cont(p); % find bif
0048 p.bifchecksw=bif; % set remaining parameters back to original
0049 p.pnamesw=pnamesw; p.ds=ds; p.nsteps=nsteps; 
0050 end % findbif
0051 
0052 % inner functions
0053 function cstop=findbiffu(p,brout,ds)
0054 cstop=0;
0055 end
0056 
0057 function findbifheadfu(p)
0058 end