星型网络中的奇异态

Chimera States in Star Networks

阅读

下载地址：

https://chinazhang-my.sharepoint.com/:b:/g/personal/rep_rebirth_zh_zhangwenhao_icu/ERB8IU5Gs8RJn2Xiea_th3cBqoBJz-l4t2Y9E0jSjN0yow?e=QnJKQK

复现

单振子

Rossler

module Rossler
    implicit none
    real,parameter :: h=0.001
    integer,parameter :: N=3,M=100,MaxT=250000,T_rans=100000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=0.284
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1
            x(k,2)=x2
            x(k,3)=x3
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: t,i,j
        real :: a,b,c,omega_i,epsilon
        a=0.15
        b=0.4
        c=8.5
        omega_i=0.41
        epsilon=0.0026
        coupling=0.0
!        do j=1,M,1
!            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
!        end do
        fx(1)=-(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(2)-xx(3)+coupling
        fx(2)=(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(1)+a*xx(2)
        fx(3)=b+xx(3)*(xx(1)-c)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网状结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    open(30,file="t_x_y_z_1.dat")
    open(40,file="t_x_y_z_2.dat")
    call x0_y0_z0()
    call star_network(M)
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_rans) then
            write(30,*) (t-T_rans)*h,x(1,1),x(1,2),x(1,3)
            write(40,*) (t-T_rans)*h,x(2,1),x(2,2),x(2,3)
        end if
    end do
    deallocate(star_matrix)
end program

• x = 0.585635185 y = -0.708050966 z = 1.42276955

Lorenz

module Rossler
    implicit none
    real,parameter :: h=0.01
    integer,parameter :: N=3,M=10,MaxT=35000,T_rans=31000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=9.0
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
!            x(k,1)=0.585635185
!            x(k,2)=-0.708050966
!            x(k,3)=1.42276955
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        coupling=0.0
!        do j=1,M,1
!            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
!        end do
        fx(1)=sigma*(xx(2)-xx(1))+coupling
        fx(2)=(r-xx(3))*xx(1)-xx(2)
        fx(3)=xx(1)*xx(2)-beta*xx(3)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    call x0_y0_z0()
    call star_network(M)
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_rans) then
            do i=1,10,1
                write(i+100,*) (t-T_rans)*h,x(i,1),x(i,2),x(i,3)
            end do
        end if
    end do
    deallocate(star_matrix)
end program

同步误差

module Rossler
    implicit none
    real,parameter :: h=0.001,PI=3.1415926
    integer,parameter :: N=3,M=2,MaxT=3500000,T_trans=1000000
    real :: x(M,N),kk=3.92
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        call random_number(x1)
        call random_number(x2)
        call random_number(x3)
        x(1,1)=0.1
        x(1,2)=0.1
        x(1,3)=0.1
        x(2,1)=0.1+0.0000001
        x(2,2)=0.1
        x(2,3)=0.1
        do k=1,M,1
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N)
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        if(i==1) then
            fx(1)=sigma*(xx(2)-xx(1))+kk*(x(2,1)-x(1,1))
            fx(2)=(r-xx(3))*xx(1)-xx(2)
            fx(3)=xx(1)*xx(2)-beta*xx(3)
        else
            fx(1)=sigma*(xx(2)-xx(1))+kk*(x(1,1)-x(2,1))
            fx(2)=(r-xx(3))*xx(1)-xx(2)
            fx(3)=xx(1)*xx(2)-beta*xx(3)
        end if

        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4



end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t
    real :: syn_error=0.0
    open(10,file="x0_y0_z0.dat")
    open(20,file="x_y_z_1.dat")
    open(30,file="x_y_z_2.dat")
    open(40,file="t_syn_error.dat")
    call x0_y0_z0()
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_trans) then
            syn_error = sqrt((x(1,1)-x(2,1))**2+(x(1,2)-x(2,2))**2+(x(1,3)-x(2,3))**2)
            write(20,*) (t-T_trans)*h,x(1,1),x(1,2),x(1,3)
            write(30,*) (t-T_trans)*h,x(2,1),x(2,2),x(2,3)
            write(40,*) (t-T_trans)*h,syn_error
        end if
    end do
end program

• 最大李指数

module Lorenz
    implicit none
    real,parameter :: h=0.001
    integer,parameter :: MaxT=2500000,N=3,M=2,T_trans=2000000
    real(kind=8) :: x(M,N),kk
contains

    subroutine fnf(xx,fx,i)
        implicit none
        real(kind=8) :: xx(N),fx(N)
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        if(i==1) then
            fx(1)=sigma*(xx(2)-xx(1))+kk*(x(2,1)-x(1,1))
            fx(2)=(r-xx(3))*xx(1)-xx(2)
            fx(3)=xx(1)*xx(2)-beta*xx(3)
        else
            fx(1)=sigma*(xx(2)-xx(1))+kk*(x(1,1)-x(2,1))
            fx(2)=(r-xx(3))*xx(1)-xx(2)
            fx(3)=xx(1)*xx(2)-beta*xx(3)
        end if

        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real(kind=8) :: x(M,N),xx(N),fx(N)
        real(kind=8) :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4


end module Lorenz



program main
    use Lorenz
    implicit none
    integer :: t,i
    real(kind=8) :: d1=0.0,d0,lam,lamda
    open(10,file="lamda.txt")
    open(20,file="d1.txt")
    d0=1.0e-7
    !初值基准点
    do i=1,1,1
        lam=0.0
        lamda=0.0
        kk=real(i)
        x(1,1)=1.0
        x(1,2)=1.0
        x(1,3)=1.0
        x(2,1)=1.0+d0
        x(2,2)=1.0
        x(2,3)=1.0
        do t=1,MaxT,1
            call rk4(x)
            d1=sqrt((x(2,1)-x(1,1))**2+(x(2,2)-x(1,2))**2+(x(2,3)-x(1,3))**2)
            write(20,*) d1
            !新的偏离点在上一次计算的两轨迹末端的连线上，且距离仍等于d0
            x(2,1)=x(1,1)+(d0/d1)*(x(2,1)-x(1,1))
            x(2,2)=x(1,2)+(d0/d1)*(x(2,2)-x(1,2))
            x(2,3)=x(1,3)+(d0/d1)*(x(2,3)-x(1,3))
            if(t>T_trans) then
                lam=lam+log(d1/d0)
            end if
        end do
        lamda=lam/((MaxT-T_trans-1)*h)
        write(10,*) kk,lamda
        write(*,*) kk,lamda
    end do
    close(10)
end program main

figure 1

module Rossler
    implicit none
    real,parameter :: h=0.01,PI=3.1415926
    integer,parameter :: N=3,M=100,MaxT=35000,T_trans=20000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=0.284
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            !            x(k,1)=0.585635185
            !            x(k,2)=-0.708050966
            !            x(k,3)=1.42276955
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: a,b,c,omega_i,epsilon
        a=0.15
        b=0.4
        c=8.5
        omega_i=0.41
        epsilon=0.0026
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
        end do
        fx(1)=-(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(2)-xx(3)+coupling
        fx(2)=(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(1)+a*xx(2)
        fx(3)=b+xx(3)*(xx(1)-c)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

    !平均相速度omega:一段时间过固定点的次数/时间段
    subroutine parameter_Omega(number,stepT)
        implicit none
        real :: omega(M)
        integer :: i,stepT,number(M)
        do i=1,M,1
            omega(i)=(2.0*PI*number(i))/real(stepT)
            write(50,*) i,omega(i)
        end do
        return
    end subroutine parameter_Omega

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,count_number(M),number(M)
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    open(50,file="omega.dat")
    count_number=0
    number=0
    call x0_y0_z0()
    call star_network(M)
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_trans) then
            do i=1,M,1
                write(i+100,*) (t-T_trans)*h,x(i,1),x(i,2),x(i,3)
                if(x(i,1)>0.0) then
                    count_number(i)=count_number(i)+1
                else if(count_number(i)>30.and.x(i,1)<0.0) then
                    number(i)=number(i)+1
                end if
            end do
        end if
    end do
    call parameter_Omega(number,MaxT-T_trans)
    deallocate(star_matrix)
end program

figure 1_A

figure 1_B

figure 2

module Rossler
    implicit none
    real,parameter :: h=0.001,PI=3.1415926
    integer,parameter :: N=3,M=100,MaxT=850000,T_trans=700000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=0.16
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            !            x(k,1)=0.585635185
            !            x(k,2)=-0.708050966
            !            x(k,3)=1.42276955
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: a,b,c,omega_i,epsilon
        a=0.15
        b=0.4
        c=8.5
        omega_i=0.41
        epsilon=0.0026
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,2)-x(i,1))
        end do
        fx(1)=-(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(2)-xx(3)+coupling
        fx(2)=(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(1)+a*xx(2)
        fx(3)=b+xx(3)*(xx(1)-c)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

    !平均相速度omega:一段时间过固定点的次数/时间段
    subroutine parameter_Omega(number,stepT)
        implicit none
        real :: omega(M)
        integer :: i,stepT,number(M)
        do i=1,M,1
            omega(i)=(2.0*PI*number(i))/real(stepT)
            write(50,*) i,omega(i)
        end do
        return
    end subroutine parameter_Omega

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,count_number(M),number(M)
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    open(50,file="omega.dat")
    count_number=0
    number=0
    call x0_y0_z0()
    call star_network(M)
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_trans) then
            do i=1,M,1
                write(i+100,*) (t-T_trans)*h,x(i,1),x(i,2),x(i,3)
                if(x(i,1)>0.0) then
                    count_number(i)=count_number(i)+1
                else if(count_number(i)>30.and.x(i,1)<0.0) then
                    number(i)=number(i)+1
                end if
            end do
        end if
    end do
    call parameter_Omega(number,MaxT-T_trans)
    deallocate(star_matrix)
end program

figure 2_A

figure 2_B

figure 3

figure 3

module Rossler
    implicit none
    real,parameter :: h=0.01
    integer,parameter :: N=3,M=100,MaxT=55000,T_trans=40000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=0.24
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: a,b,c,omega_i,epsilon
        a=0.15
        b=0.4
        c=8.5
        omega_i=0.41
        epsilon=0.0026
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,2)-x(i,1))
        end do
        fx(1)=-(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(2)-xx(3)+coupling
        fx(2)=(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(1)+a*xx(2)
        fx(3)=b+xx(3)*(xx(1)-c)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

    !同步因子R
    subroutine R_synchronization(xx,Rt)
        implicit none
        integer :: i,j,jj,t
        real :: ansRsub,x_average(MaxT),xx(M,MaxT),Rt(99,99),xxx(2,MaxT)
        do jj=2,M,1
            do j=2,M,1
                xxx=0.0
                x_average=0.0
                xxx(1,:)=xx(jj,:)
                xxx(2,:)=xx(j,:)
                do t=T_trans+1,MaxT
                    x_average(t)=sum(xxx(:,t))/real(2)
                end do
                ansRsub=0.0
                do i=1,2,1
                    ansRsub=ansRsub+sum(xxx(i,T_trans+1:MaxT)**2)/real(15000)-(sum(xxx(i,T_trans+1:MaxT))/real(15000))**2
                end do
                Rt(jj-1,j-1)=(sum(x_average(T_trans+1:MaxT)**2)/real(15000)-(sum(x_average(T_trans+1:MaxT))&
                        /real(15000))**2)/(ansRsub/real(2))
            end do
        end do

        return
    end subroutine R_synchronization

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,R_status(M-1,M-1)
    real :: xx(M,MaxT),Rt(99,99)
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    open(30,file="Rt_matrix.dat")
    call x0_y0_z0()
    call star_network(M)
    xx=0.0
    Rt=0.0
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_trans) then
            do i=1,M,1
                write(i+100,*) (t-T_trans)*h,x(i,1),x(i,2),x(i,3)
            end do
        end if
        xx(2:M,t)=x(2:M,1)
    end do
    call R_synchronization(xx,Rt)
    do i=1,M-1,1
        do j=1,M-1,1
            if(Rt(i,j)>0.99) then
                R_status(i,j)=1
            else
                R_status(i,j)=0
            end if
            write(30,*) i+1,j+1,R_status(i,j)
            write(*,*) i+1,j+1,Rt(i,j)
        end do
    end do
    deallocate(star_matrix)
end program

figure 4

module Rossler
    implicit none
    real,parameter :: h=0.01
    integer,parameter :: N=3,M=10,MaxT=35000,T_rans=31000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=9.0
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
        end do
        fx(1)=sigma*(xx(2)-xx(1))+coupling
        fx(2)=(r-xx(3))*xx(1)-xx(2)
        fx(3)=xx(1)*xx(2)-beta*xx(3)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    call x0_y0_z0()
    call star_network(M)
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_rans) then
            do i=1,10,1
                write(i+100,*) (t-T_rans)*h,x(i,1),x(i,2),x(i,3)
            end do
        end if
    end do
    deallocate(star_matrix)
end program

figure 4_A

figure 4_B

figure 5

figure 5

module Rossler
    implicit none
    real,parameter :: h=0.01
    integer,parameter :: N=3,M=100,MaxT=35000,T_rans=25000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=1.4
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,2)-x(i,1))
        end do
        fx(1)=sigma*(xx(2)-xx(1))+coupling
        fx(2)=(r-xx(3))*xx(1)-xx(2)
        fx(3)=xx(1)*xx(2)-beta*xx(3)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    open(30,file="i_t_x.dat")
    call x0_y0_z0()
    call star_network(M)
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_rans) then
            do i=1,100,1
                if(mod(t,50)==0) then
                    write(30,*) i,(t-T_rans)*h,x(i,1),x(i,2),x(i,3)
                end if    
                !write(i+100,*) (t-T_rans)*h,x(i,1),x(i,2),x(i,3)
            end do
        end if
    end do
    deallocate(star_matrix)
end program

figure 6

figure 6

module Rossler
    implicit none
    real,parameter :: h=0.001
    integer,parameter :: N=3,M=100,MaxT=260000,T_rans=220000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=5.0
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
        end do
        fx(1)=sigma*(xx(2)-xx(1))+coupling
        fx(2)=(r-xx(3))*xx(1)-xx(2)
        fx(3)=xx(1)*xx(2)-beta*xx(3)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    call x0_y0_z0()
    call star_network(M)
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_rans) then
            do i=1,M,1
                write(i+100,*) (t-T_rans)*h,x(i,1),x(i,2),x(i,3)
            end do
        end if
    end do
    deallocate(star_matrix)
end program

figure 7

figure 7

module Rossler
    implicit none
    real,parameter :: h=0.001
    integer,parameter :: N=3,M=100,MaxT=680000,T_trans=500000
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k=5.12
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
        end do
        fx(1)=sigma*(xx(2)-xx(1))+coupling
        fx(2)=(r-xx(3))*xx(1)-xx(2)
        fx(3)=xx(1)*xx(2)-beta*xx(3)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t
    real :: Syn_Error,Square_Differ,x_average,x_total
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    open(30,file="t_Syn_Error.dat")
    call x0_y0_z0()
    call star_network(M)
    do t=1,MaxT,1
        call rk4(x)
        if(t>T_trans) then
            Syn_Error=0.0
            Square_Differ=0.0
            x_average=0.0
            x_total=0.0
            do i=2,M,1
                x_total = x_total+x(i,1)
            end do
            x_average=x_total/real(M-1)
            do i=2,M,1
                Square_Differ=Square_Differ+(x(i,1)-x_average)**2
            end do
            Syn_Error=sqrt(Square_Differ/(M-1.0))
            write(*,*) (t-T_trans)*h,Syn_Error
            write(30,*) (t-T_trans)*h,Syn_Error
        end if    
    end do
    deallocate(star_matrix)
end program

figure 8

module Rossler
    implicit none
    real,parameter :: h=0.002
    integer,parameter :: N=3,M=10,MaxT=500000,T_trans=450000,cycle=10
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: a,b,c,omega_i,epsilon
        a=0.15
        b=0.4
        c=8.5
        omega_i=0.41
        epsilon=0.0026
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
        end do
        fx(1)=-(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(2)-xx(3)+coupling
        fx(2)=(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(1)+a*xx(2)
        fx(3)=b+xx(3)*(xx(1)-c)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        return
    end subroutine star_network

    !同步因子R
    subroutine R_synchronization(xx,Rt)
        implicit none
        integer :: i,j,jj,t
        real :: ansRsub,x_average(MaxT),xx(M,MaxT),Rt(M,M),xxx(2,MaxT)
        do jj=2,M,1
            do j=2,M,1
                xxx=0.0
                x_average=0.0
                xxx(1,:)=xx(jj,:)
                xxx(2,:)=xx(j,:)
                do t=T_trans+1,MaxT
                    x_average(t)=sum(xxx(:,t))/real(2)
                end do
                ansRsub=0.0
                do i=1,2,1
                    ansRsub=ansRsub+sum(xxx(i,T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(xxx(i,T_trans+1:MaxT))/&
                            real(MaxT-T_trans))**2
                end do
                if(ansRsub<0.000001) then
                    Rt(jj,j)=1.0
                else
                    Rt(jj,j)=(sum(x_average(T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(x_average(T_trans+1:MaxT))&
                            /real(MaxT-T_trans))**2)/(ansRsub/real(2))
                end if
            end do
        end do

        return
    end subroutine R_synchronization

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,ii,jj,R_status(M-1,M-1),num_syn,num_fully_syn,num_desyn,num_completely_desyn,syn_clusters,&
            num_chimera,num_syn_clusters
    real :: xx(M,MaxT),Rt(M,M)
    open(30,file="Rt_matrix.dat")
    open(40,file="k_p.dat")
    !call x0_y0_z0()
    do k=0.0,0.5,0.01
        num_fully_syn=0
        num_completely_desyn=0
        num_chimera=0
        num_syn_clusters=0
        call star_network(M)
        do jj=1,cycle,1
            call x0_y0_z0()
            xx=0.0
            Rt=0.0
            num_syn=0
            num_desyn=0
            syn_clusters=0
            do t=1,MaxT,1
                call rk4(x)
                xx(2:M,t)=x(2:M,1)
            end do
            call R_synchronization(xx,Rt)
            do i=2,M,1
                do j=2,M,1
                    write(30,*) k,jj,i,j,Rt(i,j)
                    if(i/=j) then
                        if(Rt(i,j)>0.95) then
                            num_syn=num_syn+1
                        else
                            num_desyn=num_desyn+1
                        end if
                    end if
                end do
            end do
            !同步簇
            do i=2,M,1
                do j=2,M,1
                    if(Rt(i,j)>0.99.and.j==i+1) then
                        syn_clusters=syn_clusters+1
                    end if
                end do
            end do
            !完全同步
            if(num_syn==72) then
                num_fully_syn=num_fully_syn+1
            !完全不同步
            else if(num_desyn==72) then
                num_completely_desyn=num_completely_desyn+1
            !同步簇    
            else if(syn_clusters>=5.and.num_syn/=72) then
                num_syn_clusters=num_syn_clusters+1
            !奇异态    
            else
                num_chimera=num_chimera+1
            end if
            write(*,*) jj,k,num_fully_syn,num_completely_desyn,num_chimera,num_syn_clusters
        end do
        write(40,*) k,num_fully_syn/real(cycle),num_completely_desyn/real(cycle),num_chimera/real(cycle),&
                num_syn_clusters/real(cycle)
        deallocate(star_matrix)
    end do
    close(30)
    close(40)
end program

figure 8_A

figure 8_B

figure 8_C

module Rossler
    implicit none
    real,parameter :: h=0.01
    integer,parameter :: N=3,M=100,MaxT=50000,T_trans=40000,cycle=10
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: a,b,c,omega_i,epsilon
        a=0.15
        b=0.4
        c=8.5
        omega_i=0.41
        epsilon=0.0026
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,2)-x(i,1))
        end do
        fx(1)=-(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(2)-xx(3)+coupling
        fx(2)=(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(1)+a*xx(2)
        fx(3)=b+xx(3)*(xx(1)-c)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        return
    end subroutine star_network

    !同步因子R
    subroutine R_synchronization(xx,Rt)
        implicit none
        integer :: i,j,jj,t
        real :: ansRsub,x_average(MaxT),xx(M,MaxT),Rt(M,M),xxx(2,MaxT)
        do jj=2,M,1
            do j=2,M,1
                xxx=0.0
                x_average=0.0
                xxx(1,:)=xx(jj,:)
                xxx(2,:)=xx(j,:)
                do t=T_trans+1,MaxT
                    x_average(t)=sum(xxx(:,t))/real(2)
                end do
                ansRsub=0.0
                do i=1,2,1
                    ansRsub=ansRsub+sum(xxx(i,T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(xxx(i,T_trans+1:MaxT))/&
                            real(MaxT-T_trans))**2
                end do
                if(ansRsub<0.000001) then
                    Rt(jj,j)=1.0
                else
                    Rt(jj,j)=(sum(x_average(T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(x_average(T_trans+1:MaxT))&
                            /real(MaxT-T_trans))**2)/(ansRsub/real(2))
                end if
            end do
        end do

        return
    end subroutine R_synchronization

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,ii,jj,R_status(M-1,M-1),num_syn,num_fully_syn,num_desyn,num_completely_desyn,syn_clusters,&
            num_chimera,num_syn_clusters
    real :: xx(M,MaxT),Rt(M,M)
    open(30,file="Rt_matrix.dat")
    open(40,file="k_p.dat")
    !call x0_y0_z0()
    do k=0.0,0.5,0.02
        num_fully_syn=0
        num_completely_desyn=0
        num_chimera=0
        num_syn_clusters=0
        call star_network(M)
        do jj=1,cycle,1
            call x0_y0_z0()
            xx=0.0
            Rt=0.0
            num_syn=0
            num_desyn=0
            syn_clusters=0
            do t=1,MaxT,1
                call rk4(x)
                xx(2:M,t)=x(2:M,1)
            end do
            call R_synchronization(xx,Rt)
            do i=2,M,1
                do j=2,M,1
                    write(30,*) k,jj,i,j,Rt(i,j)
                    if(i/=j) then
                        if(Rt(i,j)>0.9999) then
                            num_syn=num_syn+1
                        else
                            num_desyn=num_desyn+1
                        end if
                    end if
                end do
            end do
            !同步簇
            do i=2,M,1
                do j=2,M,1
                    if(Rt(i,j)>0.9999.and.j==i+1) then
                        syn_clusters=syn_clusters+1
                    end if
                end do
            end do
            !完全同步
            if(num_syn==9702) then
                num_fully_syn=num_fully_syn+1
            !完全不同步
            else if(num_desyn==9702) then
                num_completely_desyn=num_completely_desyn+1
            !同步簇    
            else if(syn_clusters>=10.and.num_syn/=9702) then
                num_syn_clusters=num_syn_clusters+1
            !奇异态    
            else
                num_chimera=num_chimera+1
            end if
            write(*,*) jj,k,num_fully_syn,num_completely_desyn,num_chimera,num_syn_clusters
        end do
        write(40,*) k,num_fully_syn/real(cycle),num_completely_desyn/real(cycle),num_chimera/real(cycle),&
                num_syn_clusters/real(cycle)
        deallocate(star_matrix)
    end do
    close(30)
    close(40)
end program

figure 8_D

figure 9

figure 9_A

module Rossler
    implicit none
    real,parameter :: h=0.001
    integer,parameter :: N=3,M=100,MaxT=500000,T_trans=400000,cycle=10
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            !write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,2)-x(i,1))
        end do
        fx(1)=sigma*(xx(2)-xx(1))+coupling
        fx(2)=(r-xx(3))*xx(1)-xx(2)
        fx(3)=xx(1)*xx(2)-beta*xx(3)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !        !输出矩阵
        !        do i=1,Node_number,1
        !            do j=1,Node_number,1
        !                write(20,"(F8.3)",advance='no') star_matrix(i,j)
        !            end do
        !            write(20,*)
        !        end do
        return
    end subroutine star_network

    !同步因子R
    subroutine R_synchronization(xx,Rt)
        implicit none
        integer :: i,j,jj,t
        real :: ansRsub,x_average(MaxT),xx(M,MaxT),Rt(M,M),xxx(2,MaxT)
        do jj=2,M,1
            do j=2,M,1
                xxx=0.0
                x_average=0.0
                xxx(1,:)=xx(jj,:)
                xxx(2,:)=xx(j,:)
                do t=T_trans+1,MaxT
                    x_average(t)=sum(xxx(:,t))/real(2)
                end do
                ansRsub=0.0
                do i=1,2,1
                    ansRsub=ansRsub+sum(xxx(i,T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(xxx(i,T_trans+1:MaxT))/&
                            real(MaxT-T_trans))**2
                end do
                if(ansRsub<0.000001) then
                    Rt(jj,j)=1.0
                else
                    Rt(jj,j)=(sum(x_average(T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(x_average(T_trans+1:MaxT))&
                            /real(MaxT-T_trans))**2)/(ansRsub/real(2))
                end if
            end do
        end do

        return
    end subroutine R_synchronization

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,ii,jj,R_status(M-1,M-1),num_syn,num_fully_syn,num_desyn,num_completely_desyn,syn_clusters,&
            num_chimera,num_syn_clusters
    real :: xx(M,MaxT),Rt(M,M)
    !    character(len=20) :: name1 = "t_x_y_z"
    !    character(len=10) :: name2 = ".dat"
    !    character(len=30) :: name3 = ""
    !    character(len=30) :: name(M)
    !    open(10,file="x0_y0_z0.dat")
    !    open(20,file="star_matrix.dat")
    open(30,file="Rt_matrix.dat")
    open(40,file="k_p.dat")
    !    do j=1,M,1
    !        write(name3,'(I3)') j+100
    !        write(name(j),'(A7,A3,A4)') name1,name3,name2
    !        open(j+100,file=name(j))
    !    end do
    do k=0.0,3.5,0.01
        num_fully_syn=0
        num_completely_desyn=0
        num_chimera=0
        num_syn_clusters=0
        call star_network(M)
        do jj=1,cycle,1
            call x0_y0_z0()
            xx=0.0
            Rt=0.0
            num_syn=0
            num_desyn=0
            syn_clusters=0
            do t=1,MaxT,1
                call rk4(x)
                xx(2:M,t)=x(2:M,1)
            end do
            call R_synchronization(xx,Rt)
            do i=2,M,1
                do j=2,M,1
                    write(30,*) k,jj,i,j,Rt(i,j)
                    if(i/=j) then
                        if(Rt(i,j)>0.95) then
                            num_syn=num_syn+1
                        else
                            num_desyn=num_desyn+1
                        end if
                    end if
                end do
            end do
            !同步簇
            do i=2,M,1
                do j=2,M,1
                    if(Rt(i,j)>0.99.and.j==i+1) then
                        syn_clusters=syn_clusters+1
                    end if
                end do
            end do
            !完全同步
            if(num_syn==9702) then
                num_fully_syn=num_fully_syn+1
                !完全不同步
            else if(num_desyn==9702) then
                num_completely_desyn=num_completely_desyn+1
                !同步簇
            else if(syn_clusters>=5) then
                num_syn_clusters=num_syn_clusters+1
                !奇异态
            else
                num_chimera=num_chimera+1
            end if
            write(*,*) jj,k,num_fully_syn,num_completely_desyn,num_chimera,num_syn_clusters
        end do
        write(40,*) k,num_fully_syn/real(cycle),num_completely_desyn/real(cycle),num_chimera/real(cycle),&
                num_syn_clusters/real(cycle)
        deallocate(star_matrix)
    end do
    close(30)
    close(40)
end program

figure 9_B

module Rossler
    implicit none
    real,parameter :: h=0.001
    integer,parameter :: N=3,M=100,MaxT=500000,T_trans=400000,cycle=100
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=x1*4.0-2.0
            x(k,2)=x2*4.0-2.0
            x(k,3)=x3*4.0-2.0
            !write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j
        real :: sigma,beta,r
        sigma=10.0
        beta=8.0/3.0
        r=28.0
        coupling=0.0
        do j=1,M,1
            coupling=coupling+star_matrix(i,j)*(x(j,1)-x(i,1))
        end do
        fx(1)=sigma*(xx(2)-xx(1))+coupling
        fx(2)=(r-xx(3))*xx(1)-xx(2)
        fx(3)=xx(1)*xx(2)-beta*xx(3)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
!        !输出矩阵
!        do i=1,Node_number,1
!            do j=1,Node_number,1
!                write(20,"(F8.3)",advance='no') star_matrix(i,j)
!            end do
!            write(20,*)
!        end do
        return
    end subroutine star_network

    !同步因子R
    subroutine R_synchronization(xx,Rt)
        implicit none
        integer :: i,j,jj,t
        real :: ansRsub,x_average(MaxT),xx(M,MaxT),Rt(M,M),xxx(2,MaxT)
        do jj=2,M,1
            do j=2,M,1
                xxx=0.0
                x_average=0.0
                xxx(1,:)=xx(jj,:)
                xxx(2,:)=xx(j,:)
                do t=T_trans+1,MaxT
                    x_average(t)=sum(xxx(:,t))/real(2)
                end do
                ansRsub=0.0
                do i=1,2,1
                    ansRsub=ansRsub+sum(xxx(i,T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(xxx(i,T_trans+1:MaxT))/&
                            real(MaxT-T_trans))**2
                end do
                if(ansRsub<0.000001) then
                    Rt(jj,j)=1.0
                else
                    Rt(jj,j)=(sum(x_average(T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(x_average(T_trans+1:MaxT))&
                            /real(MaxT-T_trans))**2)/(ansRsub/real(2))
                end if
            end do
        end do

        return
    end subroutine R_synchronization

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,ii,jj,R_status(M-1,M-1),num_syn,num_fully_syn,num_desyn,num_completely_desyn,syn_clusters,&
            num_chimera,num_syn_clusters
    real :: xx(M,MaxT),Rt(M,M)
!    character(len=20) :: name1 = "t_x_y_z"
!    character(len=10) :: name2 = ".dat"
!    character(len=30) :: name3 = ""
!    character(len=30) :: name(M)
!    open(10,file="x0_y0_z0.dat")
!    open(20,file="star_matrix.dat")
    open(30,file="Rt_matrix.dat")
    open(40,file="k_p.dat")
!    do j=1,M,1
!        write(name3,'(I3)') j+100
!        write(name(j),'(A7,A3,A4)') name1,name3,name2
!        open(j+100,file=name(j))
!    end do
    do k=0.0,3.5,0.01
        num_fully_syn=0
        num_completely_desyn=0
        num_chimera=0
        num_syn_clusters=0
        call star_network(M)
        do jj=1,cycle,1
            call x0_y0_z0()
            xx=0.0
            Rt=0.0
            num_syn=0
            num_desyn=0
            syn_clusters=0
            do t=1,MaxT,1
                call rk4(x)
                xx(2:M,t)=x(2:M,1)
            end do
            call R_synchronization(xx,Rt)
            do i=2,M,1
                do j=2,M,1
                    write(30,*) k,jj,i,j,Rt(i,j)
                    if(i/=j) then
                        if(Rt(i,j)>0.95) then
                            num_syn=num_syn+1
                        else
                            num_desyn=num_desyn+1
                        end if
                    end if
                end do
            end do
            !同步簇
            do i=2,M,1
                do j=2,M,1
                    if(Rt(i,j)>0.99.and.j==i+1) then
                        syn_clusters=syn_clusters+1
                    end if
                end do
            end do
            !完全同步
            if(num_syn==9702) then
                num_fully_syn=num_fully_syn+1
            end if
            !完全不同步
            if(num_desyn==9702) then
                num_completely_desyn=num_completely_desyn+1
            end if
            !奇异态
            if(num_syn>5.and.num_syn<9700) then
                num_chimera=num_chimera+1
            end if
            !同步簇
            if(syn_clusters>=5.and.num_syn/=9702) then
                num_syn_clusters=num_syn_clusters+1
            end if
            write(*,*) jj,k,num_fully_syn,num_completely_desyn,num_chimera,num_syn_clusters
        end do
        write(40,*) k,num_fully_syn/real(cycle),num_completely_desyn/real(cycle),num_chimera/real(cycle),&
                num_syn_clusters/real(cycle)
        deallocate(star_matrix)
    end do
    close(30)
    close(40)
end program

figure 10

figure 11

figure 12

figure 13

figure 14

module Rossler
    implicit none
    real,parameter :: h=0.01,PI=3.1415926
    integer,parameter :: N=3,M=10,MaxT=500000,T_trans=400000,cycle=1
    real,allocatable :: star_matrix(:,:)
    real :: x(M,N),k
contains

    subroutine x0_y0_z0()
        implicit none
        integer :: k
        real :: x1,x2,x3
        call random_seed()
        do k=1,M,1
            call random_number(x1)
            call random_number(x2)
            call random_number(x3)
            x(k,1)=4*x1-2.0
            x(k,2)=4*x2-2.0
            x(k,3)=4*x3-2.0
            write(10,*) x(k,1),x(k,2),x(k,3)
        end do
    end subroutine x0_y0_z0

    subroutine fnf(xx,fx,i)
        implicit none
        real :: xx(N),fx(N),coupling
        integer :: i,j,jj
        real :: a,b,c,omega_i,epsilon,Xm
        a=0.15
        b=0.4
        c=8.5
        omega_i=0.41
        epsilon=0.0026
        coupling=0.0
        Xm=0.0
        do j=2,M,1
            Xm=Xm+x(j,1)
        end do
        Xm=Xm/real(M-1)
        if(i==1) then
            coupling=(k/2.0)*(Xm-x(i,1))
        else
            coupling=(k/2.0)*(x(1,1)-x(i,1))
        end if
        fx(1)=-(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(2)-xx(3)+coupling
        fx(2)=(omega_i+epsilon*(xx(1)**2+xx(2)**2))*xx(1)+a*xx(2)
        fx(3)=b+xx(3)*(xx(1)-c)
        return
    end subroutine fnf


    subroutine rk4(x)
        implicit none
        integer :: i
        real :: x(M,N),xx(N),fx(N),epsilon
        real :: kx1(N),kx2(N),kx3(N),kx4(N)
        do i=1,M,1
            xx=x(i,:)
            call fnf(xx,fx,i)
            kx1=h*fx
            xx=x(i,:)+0.5*kx1
            call fnf(xx,fx,i)
            kx2=h*fx
            xx=x(i,:)+0.5*kx2
            call fnf(xx,fx,i)
            kx3=h*fx
            xx=x(i,:)+kx3
            call fnf(xx,fx,i)
            kx4=h*fx
            x(i,:)=x(i,:)+(kx1+2.0*kx2+2.0*kx3+kx4)/6.0
        end do
        return
    end subroutine rk4

    !建立网络结构
    subroutine star_network(Node_number)
        implicit none
        integer :: Node_number !节点总数
        integer :: i,j
        allocate(star_matrix(Node_number,Node_number))
        !初始化网络矩阵
        star_matrix=0
        do i=2,Node_number,1
            star_matrix(1,i)=k/2.0
            star_matrix(i,1)=k/2.0
        end do
        !输出矩阵
        do i=1,Node_number,1
            do j=1,Node_number,1
                write(20,"(F8.3)",advance='no') star_matrix(i,j)
            end do
            write(20,*)
        end do
        return
    end subroutine star_network

    !同步因子R
    subroutine R_synchronization(xx,Rt)
        implicit none
        integer :: i,j,jj,t
        real :: ansRsub,x_average(MaxT),xx(M,MaxT),Rt(M,M),xxx(2,MaxT)
        do jj=2,M,1
            do j=2,M,1
                xxx=0.0
                x_average=0.0
                xxx(1,:)=xx(jj,:)
                xxx(2,:)=xx(j,:)
                do t=T_trans+1,MaxT
                    x_average(t)=sum(xxx(:,t))/real(2)
                end do
                ansRsub=0.0
                do i=1,2,1
                    ansRsub=ansRsub+sum(xxx(i,T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(xxx(i,T_trans+1:MaxT))/&
                            real(MaxT-T_trans))**2
                end do
                if(ansRsub<0.000001) then
                    Rt(jj,j)=1.0
                else
                    Rt(jj,j)=(sum(x_average(T_trans+1:MaxT)**2)/real(MaxT-T_trans)-(sum(x_average(T_trans+1:MaxT))&
                            /real(MaxT-T_trans))**2)/(ansRsub/real(2))
                end if

            end do
        end do

        return
    end subroutine R_synchronization

    !平均相速度omega:一段时间过固定点的次数/时间段
    subroutine parameter_Omega(number,stepT)
        implicit none
        real :: omega(M)
        integer :: i,stepT,number(M)
        do i=1,M,1
            omega(i)=(2.0*PI*number(i))/real(stepT)
            write(50,*) i,omega(i)
        end do
        return
    end subroutine parameter_Omega

end module Rossler
program Rossler_Star_Networks
    use Rossler
    implicit none
    integer i,j,t,ii,jj,num_syn,num_fully_syn,num_desyn,num_completely_desyn,syn_clusters,&
            num_chimera,num_syn_clusters,count_number(M),number(M)
    real :: xx(M,MaxT),Rt(M,M)
    character(len=20) :: name1 = "t_x_y_z"
    character(len=10) :: name2 = ".dat"
    character(len=30) :: name3 = ""
    character(len=30) :: name(M)
    open(10,file="x0_y0_z0.dat")
    open(20,file="star_matrix.dat")
    open(30,file="Rt_matrix.dat")
    open(40,file="k_p.dat")
    open(50,file="omege.dat")
    do j=1,M,1
        write(name3,'(I3)') j+100
        write(name(j),'(A7,A3,A4)') name1,name3,name2
        open(j+100,file=name(j))
    end do
    do k=0.15,0.15,0.01
        num_fully_syn=0
        num_completely_desyn=0
        num_chimera=0
        num_syn_clusters=0
        call star_network(M)
        do jj=1,cycle,1
            call x0_y0_z0()
            xx=0.0
            Rt=0.0
            num_syn=0
            num_desyn=0
            syn_clusters=0
            count_number=0
            number=0
            do t=1,MaxT,1
                call rk4(x)
                do i=1,M,1
                    if(t>T_trans) then
                        write(i+100,*) (t-T_trans)*h,x(i,1),x(i,2),x(i,3)
                        if(x(i,1)>10.0) then
                            count_number(i)=count_number(i)+1
                        else if(count_number(i)>30.and.x(i,1)<10.0) then
                            number(i)=number(i)+1
                        end if
                    end if
                end do
                xx(:,t)=x(:,1)
            end do
            call parameter_Omega(number,MaxT-T_trans)
            call R_synchronization(xx,Rt)
            do i=2,M,1
                do j=2,M,1
                    write(30,*) i,j,Rt(i,j)
                    if(i/=j) then
                        if(Rt(i,j)>0.99) then
                            num_syn=num_syn+1
                        else
                            num_desyn=num_desyn+1
                        end if
                    end if
                end do
            end do
            !同步簇
            do i=2,M,1
                do j=2,M,1
                    if(Rt(i,j)>0.99.and.j==i+1) then
                        syn_clusters=syn_clusters+1
                    end if
                end do
            end do
            !完全同步
            if(num_syn==9702) then
                num_fully_syn=num_fully_syn+1
            end if
            !完全不同步
            if(num_desyn==9702) then
                num_completely_desyn=num_completely_desyn+1
            end if
            !奇异态
            if(0<num_fully_syn.and.num_fully_syn<9702) then
                num_chimera=num_chimera+1
            end if
            !同步簇
            if(syn_clusters>=2.and.num_syn/=9702) then
                num_syn_clusters=num_syn_clusters+1
            end if
        end do
        write(40,*) k,num_fully_syn/real(cycle),num_completely_desyn/real(cycle),num_chimera/real(cycle),&
                num_syn_clusters/real(cycle)
        deallocate(star_matrix)
    end do
    close(30)
    close(40)
end program

figure 14_A

figure 14_B

figure 15

figure 15_A

figure 15_B

figure 16

figure 16_A

figure 16_B

figure 16_C

figure 17