setacc.m

如何求解函数的零点的两个子函数

 %function SETACC (NZERI,NFOUND,XMIN,XMAX,IMAX,I,F,A,B)

%     ***************************************************************
%
%        L. C. ANDREANI AND A. PASQUARELLO, PISA, 1986
%
%        Bracketing of zeros: finds extrema of intervals
%        containing a specified number of zeros of function F(X).
%
%        NZERI       (ENTRY) = Number of required zeros (unmodified)
%        NFOUND      (EXIT)  = Number of zeros found
%        XMIN , XMAX (ENTRY) = extrema of bracketing region
%                              (xmin is modified if F(xmin)=0)
%        IMAX        (ENTRY) = Maximum number of steps (unmodified)
%        I           (EXIT)  = Actual number of steps
%        F=F(X)      (ENTRY) = External function (unmodified)
%        A , B       (EXIT)  = Arrays with extrema of intervals
%                              bracketing a zero
%
%     ****************************************************************
%
 %     IMPLICIT DOUBLE PRECISION (A-H,O-Z)
 %    REAL*8 A(NZERI),B(NZERI)
%
%     ****************************************************************

%             DELTA = PASSO DELLA RICERCA
%             FRACT = PARAMETRO DI AGGIUSTAMENTO
%             N = CONTATORE ZERI INDIVIDUATI
%
%     ****************************************************************
      clear all
      NZERI=30;
      %NFOUND
      XMIN=-5.;
      XMAX=5;
      IMAX=50;%能否找全零点区间跟这个值的取值有关越大越全
      %F,A,B
     FRACT=.001;
      DELTA=(XMAX-XMIN)/IMAX ;
      NFOUND=0 ;

%     CALCOLA IL PRIMO SEGNO E, SE F=0, ARRETRA

      NCONTR=0 
    if (NCONTR==4) 
         fprintf( 'FROM MODULE SETACC: FUNCTION IDENTICALLY ZERO')
         return
       end
       F1=sin(2*pi*XMIN)-.10;
      S1=sign(F1) ;
      if (F1==0.) 
         XMIN=XMIN-DELTA*FRACT ;
         NCONTR=NCONTR+1 ;
         if (NCONTR==4) 
         fprintf( 'FROM MODULE SETACC: FUNCTION IDENTICALLY ZERO')
         return
       end
      end

%     ***************************************************************
%     *  LOOP CHE SETACCIA L' INTERVALLO ASSEGNATO
%     *  CALCOLA I SEGNI AI DUE ESTREMI, E LI CONFRONTA:
%     *  SE SONO UGUALI PROSEGUE , SE SONO DIVERSI ...
%     *  ... HA TROVATO UNO ZERO, PERCIO' INCREMENTA NFOUND
%     ***************************************************************
%
      for  I=1:IMAX
         X=XMIN+I*DELTA ;
         NCONTR=0 ;
       if (NCONTR==4) 
            fprintf( 'FROM MODULE SETACC: FUNCTION IDENTICALLY ZERO')
            return
        end
         F2=sin(2*pi*X)-.10 ;
         S2=sign(F2) ;
         if (F2==0.e0) 
            NCONTR=NCONTR+1 ;
            X=X+DELTA*FRACT ;
            if (NCONTR==4) 
            fprintf( 'FROM MODULE SETACC: FUNCTION IDENTICALLY ZERO')
            return
            end
            elseif (S2~=S1) 
            NFOUND=NFOUND+1;
            A(NFOUND)=X-DELTA;
            B(NFOUND)=X ;
            if (NFOUND==NZERI)
            return
            end
            S1=S2 ;
            end
         end