pplane7.m

matlab of linear system

       set(arr.lines,'vis','on');
       set([arr.hx;arr.hy;arr.arrows;arr.barrows],'vis','off');
       set(menull,'enable','on','label','Show nullclines.');	  
      case 'arrows'
       set(arr.arrows,'vis','on');
       set([arr.lines;arr.hx;arr.hy;arr.barrows],'vis','off');
       set(menull,'enable','on','label','Show nullclines.');	  
      otherwise
       set([arr.hx;arr.hy;arr.lines;arr.arrows;arr.barrows],'vis','off');
       set(menull,'enable','on','label','Show nullclines.');	  
     end
     dud.syst.fieldtype = Arrflag;
     set(ppdisp,'user',dud);
     
     if ( (NumbFPts ~= dud.syst.npts) | (any(WINvect ~= dud.syst.wind) ) )
       dud.syst.wind = WINvect;			
       dud.syst.npts = NumbFPts;
       set(ppdisp,'user',dud);
       pplane7('dirfield',ppdisp);
     end
     figure(ppdisp);
   else
     sud.flag = 1;
     set(ppset,'user',sud);	
      
     sud = get(ppset,'user');
     %  WINvect = sud.c.wind;
     
     Xname = sud.c.xvar;
     Yname = sud.c.yvar;
     xderivstr = sud.c.xder;
     yderivstr = sud.c.yder;
     pname = sud.c.pname;
     parav = sud.c.pval;
     
     % Convert the parameters to their current values.  First remove the
     % blanks.  Also remove the periods inserted by users attempting to
     % make the function array smart.
     
     xderivstr(find(abs(xderivstr)==32))=[];
     l=length(xderivstr);
     for ( k = fliplr(findstr('.',xderivstr)))
       if (find('*/^' == xderivstr(k+1)))
	 xderivstr = [xderivstr(1:k-1), xderivstr(k+1:l)]
       end
       l=l-1;
     end
     
     yderivstr(find(abs(yderivstr)==32))=[];
     l=length(yderivstr);
     for ( k = fliplr(findstr('.',yderivstr)))
       if (find('*/^' == yderivstr(k+1)))
	 yderivstr = [yderivstr(1:k-1), yderivstr(k+1:l)];
       end
       l=l-1;
     end
     
     for kk = 1:6
       pval = parav{kk};
       if ~isempty(pval)
	 pabs = abs(pval);
	 kkk = find(pabs==32);
	 pval(kkk) = [];
	 l = length(pval);
	 for ( k = fliplr(findstr('.',pval)))
	   if (find('*/^' == pval(k+1)))
	     pval = [pval(1:k-1), pval(k+1:l)];
	   end
	   l=l-1;
	 end
	 parav{kk} = pval;
       end
     end
     
     % Build strings for the title.
     txderstr = xderivstr;
     tyderstr = yderivstr;
     
     kxder = find(abs(txderstr)==42);
     txderstr(kxder)=' '*ones(size(kxder));
     txderstr = strrep(txderstr,'-',' - ');
     txderstr = strrep(txderstr,'+',' + ');
     if (abs(txderstr(1)) == 32)
       txderstr = txderstr(2:length(txderstr));
     end
     
     kxder = find(abs(tyderstr)==42);
     tyderstr(kxder)=' '*ones(size(kxder));
     tyderstr = strrep(tyderstr,'-',' - ');
     tyderstr = strrep(tyderstr,'+',' + ');
     if (abs(tyderstr(1)) == 32)
       tyderstr = tyderstr(2:length(tyderstr));
     end
     
     tstr1 = [Xname,' '' = ', txderstr];
     tstr2 = [Yname,' '' = ', tyderstr];
     tstr = str2mat(tstr1,tstr2);
     dud.tstr = tstr;
     pstr1 = {' ';' '};
     pstr2 = {' ';' '};
     pstr3 = {' ';' '};
     
     pstring = cell(6,1);
     for kk = 1:6
       if ~isempty(parav{kk})
	 tpstr = parav{kk};
	 kxder = find(abs(tpstr)==42);  % Get rid of *s
	 tpstr(kxder)=' '*ones(size(kxder));
	 tpstr = strrep(tpstr,'-',' - ');  % Extra spaces
	 tpstr = strrep(tpstr,'+',' + ');
	 if (abs(tpstr(1)) == 32)  % Get rid of starting space
	   tpstr = tpstr(2:length(tpstr));
	 end
	 pstring{kk} = [pname{kk},' = ', tpstr];
       else
	 % pstring{kk} = [pname{kk},' = ',parav{kk}];
	 pstring{kk} = [pname{kk},' = '];
       end
     end
     
     % Get ready to do some error trapping.
     
     SS = warning;
     warning off
     XxXxXx = WINvect(1) + rand*(WINvect(2)-WINvect(1));
     YyYyYy = WINvect(3) + rand*(WINvect(4)-WINvect(3));
     err = 0;
     
     % Now we remove the backslashes (\) used to get Greek into the
     % labels. 
     txname = Xname;
     tyname = Yname;    
     xderivstr(find(abs(xderivstr)==92))=[];
     yderivstr(find(abs(yderivstr)==92))=[];
     Xname(find(abs(Xname)==92))=[];
     Yname(find(abs(Yname)==92))=[];
     
     eval([Xname,'=XxXxXx;'],'err = 1;');
     if err
       sud.flag = 0;
       set(ppset,'user',sud);
       errmsg = ['"',xstr, '" is not a valid variable name in MATLAB.'];
       fprintf('\a')
       errordlg(errmsg,'PPLANE error','on');
       return
     end
     err = 0;
     eval([Yname,'=YyYyYy;'],'err = 1;');
     if err
       sud.flag = 0;
       set(ppset,'user',sud);
       errmsg = ['"',ystr, '" is not a valid variable name in MATLAB.'];
       fprintf('\a')
       errordlg(errmsg,'PPLANE error','on');
       return
     end
     
     
     % Replace the parameters/expressions with their values.
     
     pflag = zeros(1,6);
     perr = [];
     for kk = 1:6
       if ~isempty(pname{kk})
	 pn = pname{kk};
	 pv = parav{kk};
	 if isempty(pv) 
	   perr = [perr, sud.h.pval(kk)];
	 else 
	   if isempty(str2num(pv)) % This is an expression.
	     tpv = pv;
	     tpv(find(abs(tpv)==92))=[];
	     err = 0; pval = '';
	     eval(['pval = ',tpv,';'],'err=1;');
	     if (err | isempty(pval))
	       errmsg = ['The expression for ',pn,' is not valid.'];
	       fprintf('\a')
	       errordlg(errmsg,'pplane7 error','on');
	       warning(SS)
	       return
	     end
	   end
	   xxderivstr = pplane7('paraeval',pn,pv,xderivstr);
	   yyderivstr = pplane7('paraeval',pn,pv,yderivstr);
	   if (~strcmp(xxderivstr,xderivstr) | ~strcmp(yyderivstr,yderivstr) )
	     pflag(kk) = 1;
	     xderivstr = xxderivstr;
	     yderivstr = yyderivstr;
	   end
	 end
       end
     end
     
     % We have to make the derivative strings array smart.
     
     l = length(xderivstr);
     for (k=fliplr(find((xderivstr=='^')|(xderivstr=='*')|(xderivstr=='/'))))
       xderivstr = [xderivstr(1:k-1) '.' xderivstr(k:l)];
       l = l+1;
     end
     
     l = length(yderivstr);
     for (k=fliplr(find((yderivstr=='^')|(yderivstr=='*')|(yderivstr=='/'))))
       yderivstr = [yderivstr(1:k-1) '.' yderivstr(k:l)];
       l = l+1;
     end
     
     % Some more error trapping.
     
     err = 0;res = 1;
     eval(['res = ',xderivstr, ';'],'err = 1;');
     if err | isempty(res)
        if isempty(perr)
	 errmsg = ['The first differential equation ',...
		   'is not entered correctly.'];
       else
	 errstr1 = ['The first differential equation ',...
		    'does not evaluate correctly.'];
	 errstr2 = ['At least one of the parameter values is not ',...
		    'a number.'];
	 errmsg = str2mat(errstr1,errstr2);
	 perr
	 set(perr,'string','?');
       end
       sud.flag = 0;
       set(ppset,'user',sud);
       fprintf('\a')
       errordlg(errmsg,'PPLANE error','on');
       return;
     end
     
     err = 0;res = 1;
     eval(['res = ',yderivstr, ';'],'err = 1;');
     if err | isempty(res)
       if isempty(perr)
	 errmsg = ['The second differential equation ',...
		   'is not entered correctly.'];
       else
	 errstr1 = ['The second differential equation ',...
		    'does not evaluate correctly.'];
	 errstr2 = ['At least one of the parameter values is not ',...
		    'a number.'];
	 errmsg = str2mat(errstr1,errstr2);
	 set(perr,'string','?');
       end
       sud.flag = 0;
       set(ppset,'user',sud);
       fprintf('\a')
       errordlg(errmsg,'PPLANE error','on');
       return;
     end
     
     % If an old function m-file exists delete it, and then build a new one.
     
     %    if (~strcmp(dfcn,'') & exist(dfcn)==2) delete([dfcn,'.m']);end
     tee = clock;
     tee = ceil(tee(6)*100);
     dfcn=['pptp',num2str(tee)];
     fcnstr = ['function YYyYypr = ',dfcn,'(t,YyYy)\n\n'];
     commstr = '%%%% Created by pplane7\n\n';
     varstr = [Xname,' = YyYy(1,:);', Yname,' = YyYy(2,:);\n\n'];
     lenstr = ['l = length(YyYy(1,:));\n'];
     derstr1 = ['XxXxxpr = ', xderivstr,';\n'];
     derstr2 = ['if (length(XxXxxpr) < l) XxXxxpr = XxXxxpr*ones(1,l);end\n'];
     derstr3 = ['YyYyypr = ', yderivstr,';\n'];
     derstr4 = ['if (length(YyYyypr) < l) YyYyypr = YyYyypr*ones(1,l);end\n'];
     derstr5 = 'YYyYypr = [XxXxxpr;YyYyypr];\n';
     ppf = fopen([tempdir,dfcn,'.m'],'w');
     
     fprintf(ppf,fcnstr);
     fprintf(ppf,commstr);
     fprintf(ppf,varstr);
     fprintf(ppf,lenstr);
     fprintf(ppf,derstr1);
     fprintf(ppf,derstr2);
     fprintf(ppf,derstr3);
     fprintf(ppf,derstr4);
     fprintf(ppf,derstr5);
     fclose(ppf);
     
     
     
     % Find pplane7 Display if it exists.	
     % If pplane7 Display exists, update it.  If it does not build it.
     
     ppdisp = findobj('name','pplane7 Display');
     if (~isempty(ppdisp))
       figure(ppdisp);
       dud = get(ppdisp,'user');
       dud.syst = sud.c;
       dud.settings = sud.settings;
       dfcnn = dud.function;
	 if  exist(dfcnn)==2
	   delete([tempdir,dfcnn,'.m']);
	 end
	 xmstr = [txname,' vs. t'];
	 ymstr = [tyname,' vs. t'];
	 set(dud.menu(3),'label',xmstr);
	 set(dud.menu(5),'label',ymstr);
	 menull = findobj('tag','null');
	 if ~isempty(menull)
	   delete(get(menull,'user'));
	   set(menull,'user',[]);      
	 end
     else
       ppdisp = figure('name','pplane7 Display',...
		       'numb','off',...
		       'interrupt','on',...
		       'visible','off',...
		       'tag','pplane7');
       pplane7('figdefault',ppdisp);
       dud = get(ppdisp,'user');
       dud.syst = sud.c;
       switch dud.syst.name
	case 'pendulum'
	 dud.level = 'omega^2 - 2*cos(theta)';
	otherwise
	 dud.level = ' ';
       end
       dud.settings = sud.settings;
       dud.egg = sud.egg;      
       dud.noticeflag = 1;
       dud.contours = zeros(0,1);
       fs = dud.fontsize;
       ssize = dud.ssize;
       r = ssize/10;
       ppaxw = 437*1.2;	% Default axes width
       ppaxh = 315*1.2;	% Default axes height
       ppaxl = 45*1.2;	% Default axes left
       buttw = 40*1.2;     % Default button width
       titleh = 45;    % Default title height.  This is changed later.
       nframeh = 70;   % Default notice frame height
       ppaxb = 4+nframeh+35;
       bottomedge = 38;
       ppdh = bottomedge + nframeh + ppaxh + titleh;
       uni = get(0,'units');
       set(0,'units','pixels');
       ss = get(0,'screensize');
       set(0,'units',uni);
       sw = ss(3);sh = ss(4);
       bottom = 10;
       if r*ppdh > sh -bottom -35;
	 r = (sh-bottom-35)/ppdh;
	 fs = 10*r;
	 lw = 0.5*r;
	 set(gcf,'defaultaxesfontsize',fs,'defaultaxeslinewidth',lw);
	 set(gcf,'defaulttextfontsize',fs);
	 set(gcf,'defaultlinelinewidth',lw);
	 set(gcf,'defaultuicontrolfontsize',fs*0.9);
       end
       
       % Set up the bulletin window.
       
       nframe = uicontrol('style','frame','visible','on');
       nstr = {'More';'than';'five';'lines';'of text'};
       dud.notice = uicontrol('style','text',...
			      'horiz','left',...
			      'string',nstr,'visible','on');
       ext = get(dud.notice,'extent');
       nframeh = ext(4)+2;
       
       titleh = r*titleh;
       ppaxl = r*ppaxl;
       ppaxw = r*ppaxw;
       ppaxh = r*ppaxh;
       ppaxb = nframeh+r*bottomedge;
       buttw = r*buttw;
       butth = fs+10*r;
       buttl = ppaxl + ppaxw + 5;
       buttsep = (ppaxh - butth)/2;
       % Set up the coordinate display
       
       cstr = '(0.99999999, 0,99999999)';
       dud.ccwind = uicontrol('style','text',...
			      'horiz','left',...
			      'string',cstr,...
			      'visible','on');
       cext = get(dud.ccwind,'extent');
       ccwindtxt = uicontrol('style','text',...
			     'horiz','left',...
			     'string','Cursor position: ',...
			     'visible','on');
       cwh = cext(4);
       cww = cext(3);
       
       
       % Set up the plot axes.
       
       dud.axes = axes('units','pix',...
		       'position',[ppaxl,ppaxb,ppaxw,ppaxh],...
		       'next','add',...
		       'box','on',...
		       'interrupt','on',...
		       'xgrid','on',...
		       'ygrid','on',...
		       'drawmode','fast',...
		       'visible','off',...
		       'tag','display axes');
       
       % Set up the title.
       
       dud.title.axes = axes('box','off','xlim',[0 1],'ylim',[0 1],...
			     'units','pix','vis','off',...
			     'xtick',[-1],'ytick',[-1],...
			     'xticklabel','','yticklabel','');
       
       dud.title.eq = text(0.01,0.5,' ','vert','middle');
       dud.title.p1 = text(0.75,0.5,' ','vert','middle');
       dud.title.p2 = text(0.65,0.5,' ','vert','middle');
       dud.title.p3 = text(0.55,0.5,' ','vert','middle');
       tstr = {'x_2';'y^2'};
       set(dud.title.eq,'string',tstr,'units','pix');
       ext = get(dud.title.eq ,'extent');
       titleh = ext(4)+15*r;
       set(dud.title.eq,'units','data');
       taxpos = [ppaxl,ppaxb+ppaxh,ppaxw,titleh];
       set(dud.title.axes,'pos',taxpos,'color',get(gcf,'color'));
       
       % Finish the positions.
       
       ppdw = buttl + buttw +5;
       ppdh = ppaxb+ppaxh+titleh;
       set(nframe,'pos',[10,1,ppdw-20,nframeh]);
       set(dud.notice,'pos',[15,1,ppdw-