Indeed, I have obtained the same values following the hint above. Further it would be possibly to look for bifurcations by a parameter for nonlinear and non autonomous systems? for example if I would substitute A by a function   A(t) or at least a function constant on intervals and look again for bifurcations by F1 (here g). Or more, bifurcations by two parameters at once ? Should here be a meaning?

Ok. I apologize for the mistakes along my posts. In the attach you can see my real system. Another question: after this computation, how could I do to insert the code into another file (in latex for example, but not exporting the document, I would like with copy/paste) just as it is.  More precisely i want to appear 0.0012 not something times with "e", the Jacobian as a matrix and so on. Thanks a lot!

p.s I have Maple 9.5 and start a sheet running Maple 9.5 not Classic Worksheet 9.5. I note that it runs rather slowly.

View 9134_bifurcAcaseII.mw on MapleNet or Download 9134_bifurcAcaseII.mw
View file details

Sorry, pars:= {b =1, c=0.2, d=0.5, e=14, f=0.15, g=2.2} (my original system has coefficients like E2, E2..etc) I studied the system using the commands you proposed above; it worked for e, g but as I see not for a.

I'm back. 1. Trying to look for bifurcations by a, for the set of values pars:= {F2 =1, E2=0.2, D2=0.5, D1=14, E1=0.15, F1=2.2} Using the command > evalf([allvalues(eval(S,pars))]) I get the the following Error message, (in evala/preproc3) floats not handled yet. What's wrong? 2. Is it posiibly to find two bifurcations points for the same parameter?

It worked! Some questions 1. Where do appear the initial conditions x(0)=y(0)=0.01? 2. the bifurcation is considered at e=85.4 or at e=85.5? Where should I stop, where are still two equillibria points or at least one or no more one? Well, I tried for my real system (actually I thought is a slight modification but now I think it is significantly modified) So x'=a*b*((x^2+c^2)/(x^2+1))-d*x y'=a*e*((y^3+f^3)/(y^3+1)*(g*x+1))-y (a modification appears in the down term, for the first bracket I have y instead of x ) I followed the above approach with the following changes pars:= { a=7, b = 0.2, c=0.2, d=0.1, f=0.5, g=2.8}; and I got >evalf([allvalues(eval(S,pars))]); [{e = -0.1155507917e-1-.2044789251*I, x = 0.3444817452e-1-.1975123882*I, y = -.5608636000-.9714442513*I}, {e = -.1713064041-.1122464547*I, x = 0.3444817452e-1+.1975123882*I, y = -.5608636000-.9714442513*I}, {e = -0.3143281871e-2-0.5444323904e-2*I, x = 13.93110365, y = -.5608636000-.9714442513*I}, {e = -.1713064041+.1122464547*I, y = -.5608636000+.9714442513*I, x = 0.3444817452e-1-.1975123882*I}, {e = -0.1155507917e-1+.2044789251*I, x = 0.3444817452e-1+.1975123882*I, y = -.5608636000+.9714442513*I}, {x = 13.93110365, y = -.5608636000+.9714442513*I, e = -0.3143281871e-2+0.5444323904e-2*I}, {e= .1828614833+0.9223247044e-1*I, y = 1.121727200, x = 0.3444817452e-1-.1975123882*I}, {e = .1828614833-0.9223247044e-1*I, x = 0.3444817452e-1+.1975123882*I, y = 1.121727200}, {e = 0.006286563743, x = 13.93110365, y = 1.121727200}, {e = -0.1491201246e-1-.2638832880*I, y = -.2228705874-.3860231810*I, x =.3444817452e-1-.1975123882*I}, {x = 0.3444817452e-1+.1975123882*I, e = -.2210736249-.1448558256*I, y = -.2228705874-.3860231810*I}, {x = 13.93110365, y = -.2228705874-.3860231810*I, e = -0.4056454986e-2-0.7025986134e-2*I}, {e = -.2210736249+.1448558256*I, y = -.2228705874+.3860231810*I, x = .3444817452e-1-.1975123882*I}, {x = 0.3444817452e-1+.1975123882*I, e = -0.1491201246e-1+.2638832880*I, y = .2228705874+.3860231810*I}, {x = 13.93110365, y = -.2228705874+.3860231810*I, e = -0.4056454986e-2+0.7025986134e-2*I}, {e = .2359856373+.1190274624*I, y = .4457411749, x = 0.3444817452e-1-.1975123882*I}, {e= .2359856373-.1190274624*I, x = 0.3444817452e-1+.1975123882*I, y = .4457411749}, {e = 0.008112909971, x = 13.93110365, y = .4457411749}] In this case I have two possibilities: {e = 0.006286563743, x = 13.93110365, y = 1.121727200}, {e = 0.008112909971, x = 13.93110365, y = .4457411749}] for e=0.006 I obtain no equilibria points. With display([implicitplot(e006,x=13.92 .. 13.96,y=1.1..1.3,colour=[green,blue]), DEplot(zip(`=`,[D(x)(t),D(y)(t)],eval(D006,{x=x(t),y=y(t)})),[x(t),y(t)],t=1..20,x=13.92..13.96, y=1.1..1.3,arrows=slim)]); I can see the green line; If I take x=13.92..13.95 no line I can see, how do I interprete this fact? for e=0.008 again I do not obtain any equilibria point, but for e=0.0081, I do; what is the conclusion? there is no bifurcation? :(( Where I'm wrong? Thanks for your patience!

Robert, there were of very help for me your explanations, but I have some questions. I'm a beginner in Maple, so I'm reading all commands step by step (working with Maple 9.5). I got an error massage after the command > S:= solve([op(F), J[2,2]], [x,y,e]); I understand that you compute J[2,2] (a proper value) at the equilibrium point. OK, I deduce I should try the same algorithm for J[1,1], isn't it? Now remain only two questions: why doesn't work the above command?maybe because working with Maple 9.5 :( and what means _Z? Thanks again you all! carmina