%
% << noeudcol >> 
%
function noeudcol
clear all
close all
global mu alpha seuil
%
disp('Bifurcation noeud-col. Patientez un peu SVP !')
% parametres 
tmax=3; nt=10;  nx0=11; nmu=11;
%tmax=1; nt=10;  nx0=4; nmu=4;
%tmax=1; nt=10;  nx0=4; nmu=11;
%
alpha=1;
tmin=0;
dt=(tmax-tmin)/(nt-1);
tspan=[tmin:dt:tmax];
%%%
% initialize x
x0max=2; x0min=-2;
dx0=(x0max-x0min)/(nx0-1);
x0v=x0min:dx0:x0max;
% seuil abs(x) < seuil pour calcul f(x)
seuil=x0max*1.5;
% initialize mu
mumax=1; mumin=-.2;
dmu=(mumax-mumin)/(nmu-1);
muv=mumin:dmu:mumax;
%
% initialize mup
mup=0:.01:mumax;
%
% figure d'attente
plot(mup,sqrt(mup),'r','LineWidth',2)
axis([mumin-.1 mumax+.1 x0min-.1 x0max+.1])
grid
hold on
plot(mup,-sqrt(mup),'--r','LineWidth',2)
hold on
ylabel('x','FontSize',20)
xlabel('\mu','FontSize',20)
for imu=1:nmu
	mu=muv(imu);
	for ix0=1:nx0
		plot(muv(imu),x0v(ix0),'o','Color','k','LineWidth',2)
	end
end
disp('Changement de taille de l''image pendant l''attente !')
drawnow
%pause
%% Resolution ode
nla=0; nlamax=nt-1;
for imu=1:nmu
	for ix0=1:nx0
		mu=muv(imu);
		x0=x0v(ix0);
		[t,x]=ode23('nc',tspan,x0);
		for it=1:nt
			xa(imu,ix0,it)=x(it);
		end
	end
end
%% Boucle temporelle
hold off
for it=1:nt
	plot(mup,sqrt(mup),'r','LineWidth',2)
	axis([mumin-.1 mumax+.1 x0min-.1 x0max+.1])
	grid
	hold on
	plot(mup,-sqrt(mup),'--r','LineWidth',2)
	hold on
	ylabel('x','FontSize',20)
	xlabel('\mu','FontSize',20)
	if it==1
		for imu=1:nmu
			for ix0=1:nx0
				plot(muv(imu), xa(imu,ix0,it),'o','Color','k','LineWidth',2)
			end
		end
	elseif  it==2
		for imu=1:nmu
			for ix0=1:nx0
				muvimu=muv(imu); mudeux=[muvimu muvimu];
				xadeux=[xa(imu,ix0,it) xa(imu,ix0,it-1)];
				plot(mudeux, xadeux,'-b')
				plot(muv(imu), xa(imu,ix0,it),'o','Color','k','LineWidth',2)
				plot(muv(imu), xa(imu,ix0,it-1),'o','Color','b','LineWidth',1)
			end
		end
	else
		for imu=1:nmu
			for ix0=1:nx0
				muvimu=muv(imu); mutrois=[muvimu muvimu muvimu];
				xatrois=[xa(imu,ix0,it) xa(imu,ix0,it-1) xa(imu,ix0,it-2)];
				plot(mutrois, xatrois,'-b')
				plot(muv(imu), xa(imu,ix0,it),'o','Color','k','LineWidth',2)
				plot(muv(imu), xa(imu,ix0,it-1),'o','Color','b','LineWidth',1)
				plot(muv(imu), xa(imu,ix0,it-2),'o','Color','b','LineWidth',1)
			end
		end
	end
	pause(0)
	hold off
	if nla==0 
		% disp(' You can change the window size now !')
		% pause
		Mla=moviein(nlamax);
	end
	nla=nla+1;
	Mla(:,nla)=getframe; 
	%filew(gcf,it);
end
%%% N times, FPS frame per second
%pause
movie(Mla,5,5)
return
