Sources Esope | Cast3M

laguer
C LAGUER    SOURCE    PV        07/11/23    21:17:36     5978
      subroutine laguer(a,m,x,its,kerre)
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
      dimension frac(8)
      real*8 err,abp,abm,abx,frac
      complex dx,x1,b,d,f,g,h,sq,gp,gm,g2
      complex a(m+1),x
      data frac/.5,.25,.75,.13,.38,.62,.88,1./
*      save frac
*      write(6,*) ' entree dans laguer m ',m
      epss = 2.e-7
      mr = 8
      mt=10
      maxit = mt*mr*3
      do 1 iter = 1,maxit
      its=iter
      b=a(m+1)
      err = abs(b)
      d=cmplx(0.,0.)
      f=cmplx(0.,0.)
      abx = abs(x)
      do 2 j=m,1,-1
          f = x*f+d
          d=x*d+b
          b=x*b+a(j)
          err=abs(b)+abx*err
    2  continue
       err=epss*err
      if(abs(b).le.err) then
         return
      else
         g=d/b
         g2=g*g
         h=g2-2.*f/b
         sq=sqrt((m-1)*(m*h-g2))
         gp=g+sq
         gm=g-sq
         abp=abs(gp)
         abm=abs(gm)
         if(abp.lt.abm) gp=gm
         if(max(abp,abm).gt.0.) then
           dx=m/gp
         else
           dx=exp(cmplx(log(1.+abx),float(iter)))
         endif
      endif
      x1=x-dx
      if(x.eq.x1) return
      if(mod(iter,MT).ne.0) then
        x=x1
      else
        x=x-dx*frac(iter/mt)
      endif
   1  continue
      kerre =944
      return
      end