非线性振动奇点的求解及其类型的判断

博大广阔

%求解系统的极点和并判断其类型
%求解的系统类型：
% dx1=f(x1,x2);
% dx2=g(x1,x2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function CC=critical_points()
clc
clear all
syms x1 x2
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%5
[F,V Cri]=aa();                                               %奇点  
%F是两个函数组成的矩阵，V是函数中得变量V=[x1 x2];Cri为奇点
J=jacobin(F,V);
for k=1:1:size(Cri,1)
A=subs(J,{x1,x2},{Cri(k,1),Cri(k,2)});
A=double(A);
p=trace(A);
q=det(A) ;       %节点类型取决于A的特征值结构
w=p^2-4*q;
%奇点类型判断%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if w>=0                                    %w>0,r1和r2 是不等实根，w=0是重根
    if q>0                                   %R1和R2是同号
        if p<=0
            Cr=[Cri(k,1),Cri(k,2)]
            disp('稳定节点')
        elseif p>0
            Cr=[Cri(k,1),Cri(k,2)]
            disp('不稳定节点')
        end
            
    elseif q<0
         Cr=[Cri(k,1),Cri(k,2)]
         disp('鞍点')
    end
      
elseif w<0                  %共轭复根
  if p==0  
      Cr=[Cri(k,1),Cri(k,2)]
      disp('中点')
  else p~=0
      if p<=0
         Cr=[Cri(k,1),Cri(k,2)]
          disp('稳定焦点')         
      elseif p>0
          Cr=[Cri(k,1),Cri(k,2)]
          disp('不稳定焦点')            
      end
  end
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
end
end
%%%%%%%%%%%不同的函数需在此设置
function [F,V Cri]=aa()
%F是两个函数组成的矩阵，V是函数中得变量V=[x1 x2];Cri为奇点
syms x1 x2
F=[x1*(1-x1-x2);0.25*x2*(2-3*x1-x2)];%两个函数就是两行
V=[x1,x2];  
[x1,x2]=solve('x1*(1-x1-x2)=0','0.25*x2*(2-3*x1-x2)= 0'); %函数，不同的需更改
Cri=[x1,x2];
end