      subroutine nyra(n,u,up,t,tend,rtol,atol,f2nd)
c
c-----!-----------------------------------------------------------------
c
c   NYRA ... A Runge-Kutta-Nystroem-type block predictor-
c            corrector method for parallel computers with shared memory
c
c
c     NYRA solves the second order initial value problem for nonstiff
c     ODEs of the form
c 
c        2
c     d u
c     ---  = f2nd (t, u),          u=(u(1),..,u(n)
c        2
c     d t
c 				    1  	n     n
c 			     f2nd: R x R --> R
c 
c 
c     with a 7-order variable stepsize 4-stage method based on direct
c     collocation at Radau-IIA points.  The implicit relations in the
c     corrector steps are solved by parallel fix-point iterations.
c       
c     For the predictor, the code computes approximations with four 
c     different stepsizes, hence the maximal number of processors is 
c     4 x 4 = 16.
c 
c     Version of September 22, 1998
c     ( This is an experimental version. )
c 
c     The code contains special directives for the Convex-SPP
c     Fortran-Compiler of the form "c$dir loop_parallel".  These
c     statements are comments for other compilers and will be ignored,
c     the code runs in sequential mode.
c 
c
c  Authors:
c            
c     N. H. Cong,         
c     Faculty of Mathematics, Mechanics and Informatics 
c     Hanoi University of Sciences
c     90 Nguyen Trai, Dong Da, Hanoi, Vietnam
c     e-mail: nhcong@ncst.ac.vn
c           
c     H. Podhaisky, R. Weiner, K. Strehmel
c     Fachbereich Mathematik und Informatik, Universitaet Halle
c     06099 Halle, Germany
c     e-mail: podhaisky@mathematik.uni-halle.de
c             rw@mail.mathematik.uni-halle.de
c	      strehmel@.mathematik.uni-halle.de
c
c
c Please do not hesitate to send us remarks, errors and suggestions
c by e-mail to: 
c	    podhaisky@mathematik.uni-halle.de 	
c
c-----!-----------------------------------------------------------------
c    
c  References:
c
c	[1] N.H. Cong, K. Strehmel and R. Weiner:
c	Runge-Kutta-Nystroem-type parallel block predictor-corrector
c	methods. Tech. Report 29 (1997), University Halle
c
c
c For retrieving the latest version visit URL:
c http://www.mathematik.uni-halle.de/institute/numerik/software
c
c-----!----------------------------------------------------------------- 
c
c Parameters:
c 
c  n    =  Number of equations of the ODE system
c
c  f2nd = Name (external) of the user-supplied subroutine computing
c	  the right hand side of the ODE. It must have the form
c
c	  subroutine f2nd(n,t,u,fw)
c         integer n
c         real*8 t,u(n),fw(n)
c    		
c		...
c
c  u   =   On Input:  initial values for u.
c          On Output: solution u in t=tend.
c
c  up  =   On Input:  initial values for the first derivative of u
c          On Output: derivative in t=tend
c
c
c   t  =   The time-variable "t" in the ODE. On exit, t should be tend.
c
c tend =   Final point for integration.
c
c
c nx   =   A compiling time constant (see code below). Be careful:
c	   nx must be greater or equal n !
c
c lpirkn = parameter in common block below:
c                  logical lpirkn
c                  common /clpirkn /lpirkn
c
c          .true.  -- use PIRKN method (without block predictor)
c          .false. -- use BIKRN method (block predictor)
c
c-----!----------------------------------------------------------------- 

      implicit none

      integer nx
      parameter (nx=10000)
      integer n

      integer s,r
      parameter(s=4,r=s)
      real*8 an(s,s),c(s),aa(s),bn(s),b(s),vss(s,s,s),wss(s,s,s)
      real*8 wss0(s,s,r),vss0(s,s,r)
      real*8 be(s),bne(s)
      real*8 fuu(nx,s,r)

      real*8 tmp(s,2*s), y(s,2*s), xinv(2*s,2*s)

      real*8 u(n),up(n),t,tend,rtol,atol,h,hs,h2,hnew
      real*8 un(nx,s,2),uround,err,uex(nx)
      real*8 scal(nx),ue(nx,r)

      real*8 dnrm2,hold

      integer icur,inext

      real*8 upn(nx,r,2)

      real*8 alpha

      data uround /1d-15/
      external f2nd
      logical fixed
      data fixed /.false./
      integer nfeval,nsteps,nrejct
      common /cstat/ nsteps,nrejct,nfeval

      integer nchange
      common /cnchange/nchange

      integer ns,i


      call coeff (s,c,an,bn,b,aa,vss0,wss0,be,bne)

      alpha=1d0

      icur=1
      inext=3-icur

      do i=1,r
         call dcopy(n,u,1,un(1,i,icur),1)
         call dcopy(n,up,1,upn(1,i,icur),1)
      end do

      nsteps=0
      nfeval=0
      nrejct=0
      nchange=0
c
c choose a fixed inital stepsize h
c
      h=(tend-t)/50d0/(-1d0*dlog10(1d0*atol))


c
c begin integration
c

 1000 continue
      nsteps=nsteps+1
      if (nsteps.gt.50000) then
         write (*,*) 'Too many steps'
         return
      end if
      do  i=1,n
         scal(i)=1d0/(atol+rtol*dabs(un(i,1,icur)))
      end do

      if (nsteps-nrejct.gt.1) then
         alpha=h/hold
      end if

      if (abs(alpha-1d0).gt.1d-5) then
         nchange=nchange+1
         call recoeff (s,aa,c,alpha,vss,wss,tmp,y,xinv)

         call step (n,t,h,un(1,1,icur),upn(1,1,icur),
     $        un(1,1,inext),upn(1,1,inext),ue,
     $        c,an,bn,be,bne,b,aa,vss,wss,fuu,alpha)

      else

         call step (n,t,h,un(1,1,icur),upn(1,1,icur),
     $        un(1,1,inext),upn(1,1,inext),ue,
     $        c,an,bn,be,bne,b,aa,vss0,wss0,fuu,alpha)

      end if

      if (.not.fixed) then
         err=0d0
         call daxpy(n,-1d0,un(1,1,inext),1,ue(1,1),1)
         do i =1,n
            err=err+(ue(i,1)*scal(i))**2
         end do
         err=dsqrt(err/n)
         hnew=dmin1(3d0,dmax1(0.3d0,(1d0/err)**(1d0/(s+1))*0.8))
         if (.true.) then
            if (0.90d0.lt.hnew.and.hnew.lt.2d0) then
               hnew=1d0
               err=0.9d0
            end if
         end if
         hnew=hnew*h
      else
         err=0.9d0
         hnew=h
      end if

      if (err.lt.1d0 )then
         icur=3-icur
         inext=3-inext
         t=t+h
      else
         nrejct=nrejct+1
      end if
      if (t.lt.tend-1d3*uround) then
         hnew=min(tend-t,hnew)
         hold=h
         h=hnew
         if (hnew.lt.1d-12) then
            write (*,*) 'stepsize too small!',h
            return
         end if
         goto 1000
      end if

      call dcopy(n,un(1,1,icur),1,u,1)

      return
      end


c-----!-----------------------------------------------------------------
c --- S U B R O U T I N E   S T E P
c-----!-----------------------------------------------------------------
c
c
      subroutine step (n,t,h,u,up,un,upn,ue,c,
     $     an,bn,be,bne,b,aa,vss,wss,fuu,alpha)
      implicit none

      integer n,s,nx,r
      parameter (nx=10000,s=4,r=s)

      real*8 t,u(nx,r),up(nx,r),un(nx,r),upn(nx,r),h
      real*8 ue(nx,r)
      real*8 c(s),an(s,s),bn(s),be(s),bne(s),b(s),aa(r)
      real*8 vss(s,s,r),wss(s,s,r)

      real*8 err,uex(nx),alpha
      real*8 dnrm2

      real*8 uu(nx,s,r,2), fuu(nx,s,r)
      integer icur,inext,m

      parameter (m=1)

      integer i,j,l,k,iter,irr

      integer nfeval,nsteps,nrejct
      common /cstat/ nsteps,nrejct,nfeval

      logical lpirkn
      common /clpirkn /lpirkn

      integer iim,nths,ith


      icur=1
      inext=3-icur

      irr=r
      if (lpirkn) then
         irr=1
      end if

      iim=m
      if (nsteps-nrejct.eq.0) iim=max(m,5)
      if (abs(alpha-1d0).gt.1d-5) iim=max(m,3)
      nfeval=nfeval+irr*s*(1+iim)


      nths=irr*s

c$dir loop_private(k,i,j)
c$dir loop_parallel
      do ith=1,nths
         k=(ith+s-1)/r
         i=mod(ith,r)+1

c     Predictor

         if (nsteps-1.eq.nrejct.or.lpirkn) then
            call dcopy(n,u(1,1),1,uu(1,i,k,icur),1)
            call daxpy(n,aa(k)*h*c(i),up(1,1),1,uu(1,i,k,icur),1)
         else
            call dcopy(n,0d0,0,uu(1,i,k,icur),1)
            do j=1,s
               call daxpy(n,vss(i,j,k),u(1,j),1,uu(1,i,k,icur),1)
               call daxpy(n,h*wss(i,j,k),up(1,j
     +              ),1,uu(1,i,k,icur),1)
            end do
         end if
         call f2nd(n,t+aa(k)*c(i)*h,uu(1,i,k,icur),fuu(1,i,k))
      end do

c pirk iteration

      do iter=1,iim
c$dir loop_private(k,i,j)
c$dir loop_parallel
         do ith=1,nths
            k=(ith+s-1)/r
            i=mod(ith,r)+1

            call dcopy(n,u(1,1),1,uu(1,i,k,inext),1)
            call daxpy(n,aa(k)*h*c(i),up(1,1),1,uu(1,i,k,inext),1)
            do j=1,s
               call daxpy(n,(aa(k)*h)**2*an(i,j),
     $              fuu(1,j,k),1,uu(1,i,k,inext),1)
            end do
         end do

         icur=(-1*(-1)**iter+3)/2
         inext=3-icur
c$dir loop_private(k,i)
c$dir loop_parallel
         do ith=1,nths
            k=(ith+s-1)/r
            i=mod(ith,r)+1
            call f2nd(n,t+aa(k)*c(i)*h,
     +           uu(1,i,k,icur),fuu(1,i,k))
         end do
      end do
c
c compute new u
c

c$dir loop_private(i)
c$dir loop_parallel
      do k=1,irr
         call dcopy(n,u(1,1),1,un(1,k),1)
         call daxpy(n,aa(k)*h,up(1,1),1,un(1,k),1)
         do i=1,s
            call daxpy(n,(aa(k)*h)**2*bn(i),fuu(1,i,k),1,un(1,k),1)
         end do
c
c compute new u'
c

         call dcopy(n,up(1,1),1,upn(1,k),1)
         do i=1,s
            call daxpy(n,aa(k)*h*b(i),fuu(1,i,k),1,upn(1,k),1)
         end do
      end do

c
c compute embedded solution
c
      call dcopy(n,u(1,1),1,ue(1,1),1)
      call daxpy(n,h,up(1,1),1,ue(1,1),1)
      do i=1,s
         call daxpy(n,h*h*bne(i),fuu(1,i,1),1,ue(1,1),1)
      end do

      return
      end

c
c-----!-----------------------------------------------------------------
c --- S U B R O U T I N E   C O E F F
c-----!-----------------------------------------------------------------
c
c
      subroutine coeff(s,c,an,bn,b,aa,vss,wss,be,bne)
      implicit none

      integer s
      real*8 c(s),an(s,s),bn(s),b(s),aa(s),vss(s,s,s),wss(s,s,s)
      real*8 be(s),bne(s)


      goto (100,200,300,400,500,600) s
      write (*,*) 'illegal s=',s
      stop "coeff"

 100  continue
      stop "s=1 in coeff"

 200  continue

      c(1) = 0.3333333333333333D0
      c(2) = 0.1D1
      an(1,1) = 0.1111111111111111D0
      an(1,2) = -0.5555555555555556D-1
      an(2,1) = 0.5D0
      an(2,2) = 0.1D-99
      bn(1) = 0.5D0
      bn(2) = 0.1D-99
      b(1) = 0.75D0
      b(2) = 0.25D0
      aa(1) = 1.D0
      aa(2) = 0.2D1
      be(1) = 0.D0
      be(2) = 0.1D1
      bne(1) = 0.D0
      bne(2) = 0.1D1

      vss(1,1,1) = 0.7407407407407407D0
      vss(1,1,2) = 0.2592592592592593D0
      vss(1,2,1) = 0.2592592592592593D0
      vss(1,2,2) = 0.7407407407407407D0
      vss(2,1,1) = 0.D0
      vss(2,1,2) = 0.5D1
      vss(2,2,1) = 0.1D1
      vss(2,2,2) = -0.4D1
      wss(1,1,1) = 0.1481481481481481D0
      wss(1,1,2) = 0.7407407407407407D-1
      wss(1,2,1) = -0.7407407407407407D-1
      wss(1,2,2) = -0.1481481481481481D0
      wss(2,1,1) = 0.D0
      wss(2,1,2) = 0.2D1
      wss(2,2,1) = 0.D0
      wss(2,2,2) = 0.4D1

      return

 300  continue

      an(1,1) = 0.1637627564304205D-1
      an(1,2) = -0.6866522903509781D-2
      an(1,3) = 0.2510657549140603D-2
      an(2,1) = 0.1811998562368431D0
      an(2,2) = 0.2862372435695795D-1
      an(2,3) = -0.1843990882473936D-2
      an(3,1) = 0.3180413817439772D0
      an(3,2) = 0.1819586182560228D0
      an(3,3) = 0.1D-99
      bn(1) = 0.3180413817439772D0
      bn(2) = 0.1819586182560228D0
      bn(3) = 0.1666666666666667D-100
      c(1) = 0.1550510257216822D0
      c(2) = 0.6449489742783178D0
      c(3) = 1.D0
      b(1) = 0.3764030627004673D0
      b(2) = 0.5124858261884216D0
      b(3) = 0.1111111111111111D0

      bne(1) = 0.D0
      bne(2) = 0.9388321936425753D0
      bne(3) = -0.4388321936425753D0

      aa(1) = 1.D0
      aa(2) = 1d0+c(2)
      aa(3) = 1d0+c(3)

      vss(1,1,1) = 1.D0
      vss(1,1,2) = 1.D0
      vss(1,1,3) = 1.D0
      vss(1,2,1) = 0.D0
      vss(1,2,2) = 0.D0
      vss(1,2,3) = 0.D0
      vss(1,3,1) = 0.D0
      vss(1,3,2) = 0.D0
      vss(1,3,3) = 0.D0
      vss(2,1,1) = 1.D0
      vss(2,1,2) = 1.D0
      vss(2,1,3) = 1.D0
      vss(2,2,1) = 0.D0
      vss(2,2,2) = 0.D0
      vss(2,2,3) = 0.D0
      vss(2,3,1) = 0.D0
      vss(2,3,2) = 0.D0
      vss(2,3,3) = 0.D0
      vss(3,1,1) = 1.D0
      vss(3,1,2) = 1.D0
      vss(3,1,3) = 1.D0
      vss(3,2,1) = 0.D0
      vss(3,2,2) = 0.D0
      vss(3,2,3) = 0.D0
      vss(3,3,1) = 0.D0
      vss(3,3,2) = 0.D0
      vss(3,3,3) = 0.D0
      wss(1,1,1) = 0.1550510257216822D0
      wss(1,1,2) = 0.2550510257216822D0
      wss(1,1,3) = 0.3101020514433644D0
      wss(1,2,1) = 0.D0
      wss(1,2,2) = 0.D0
      wss(1,2,3) = 0.D0
      wss(1,3,1) = 0.D0
      wss(1,3,2) = 0.D0
      wss(1,3,3) = 0.D0
      wss(2,1,1) = 0.6449489742783178D0
      wss(2,1,2) = 0.1060908153700972D1
      wss(2,1,3) = 0.1289897948556636D1
      wss(2,2,1) = 0.D0
      wss(2,2,2) = 0.D0
      wss(2,2,3) = 0.D0
      wss(2,3,1) = 0.D0
      wss(2,3,2) = 0.D0
      wss(2,3,3) = 0.D0
      wss(3,1,1) = 1.D0
      wss(3,1,2) = 0.1644948974278318D1
      wss(3,1,3) = 2.D0
      wss(3,2,1) = 0.D0
      wss(3,2,2) = 0.D0
      wss(3,2,3) = 0.D0
      wss(3,3,1) = 0.D0
      wss(3,3,2) = 0.D0
      wss(3,3,3) = 0.D0

      return

 400  continue

c
c indirect collocation
c

      c(1) = 0.8858795951270395D-1
      c(2) = 0.4094668644407347D0
      c(3) = 0.7876594617608471D0
      c(4) = 0.1D1
      an(1,1) = 0.6728344411839956D-2
      an(1,2) = -0.6260584380138599D-2
      an(1,3) = 0.6465233728895639D-2
      an(1,4) = -0.3009080475284759D-2
      an(2,1) = 0.681456600536828D-1
      an(2,2) = 0.2015031261381731D-1
      an(2,3) = -0.7623261099229043D-2
      an(2,4) = 0.315884496919244D-2
      an(3,1) = 0.1553032519913779D0
      an(3,2) = 0.1426744060777562D0
      an(3,3) = 0.1393766950495498D-1
      an(3,4) = -0.1711613723395395D-2
      an(4,1) = 0.2009319137389596D0
      an(4,2) = 0.2292411063595862D0
      an(4,3) = 0.6982697990145412D-1
      an(4,4) = 0.421D-29
      bn(1) = 0.2009319137389596D0
      bn(2) = 0.2292411063595862D0
      bn(3) = 0.6982697990145412D-1
      bn(4) = 0.421D-29
      b(1) = 0.2204622111767684D0
      b(2) = 0.3881934688431719D0
      b(3) = 0.3288443199800597D0
      b(4) = 0.625D-1

      be(1) = 0.D0
      be(2) = 0.1017139272442552D1
      be(3) = -0.4740236823861008D0
      be(4) = 0.456884409943549D0

      bne(1) = 0.D0
      bne(2) = 0.8697416981467024D0
      bne(3) = -0.84900396455074D0
      bne(4) = 0.4792622664040376D0

c
c direct collocation
c

      if (.true.) then
      an(1,1) = 0.5382675529471965D-2
      an(1,2) = -0.2421591783257662D-2
      an(1,3) = 0.1564645634154885D-2
      an(1,4) = -0.6018160950569517D-3
      an(2,1) = 0.6955830402056265D-1
      an(2,2) = 0.1612025009105385D-1
      an(2,3) = -0.2478766567991471D-2
      an(2,4) = 0.631768993838488D-3
      an(3,1) = 0.1545378137303791D0
      an(3,2) = 0.1448580872610296D0
      an(3,3) = 0.1115013560396398D-1
      an(3,4) = -0.3423227446790789D-3
      an(4,1) = 0.2009319137389596D0
      an(4,2) = 0.2292411063595862D0
      an(4,3) = 0.6982697990145412D-1
      an(4,4) = 0.1D-99


      bn(1) = 0.2009319137389596D0
      bn(2) = 0.2292411063595862D0
      bn(3) = 0.6982697990145412D-1
      bn(4) = 0.2333333333333333D-99


      bne(1) = 0.D0
      bne(2) = 0.8024698845740193D0
      bne(3) = -0.6619165843770475D0
      bne(4) = 0.3594466998030282D0


      c(1) = 0.8858795951270395D-1
      c(2) = 0.4094668644407347D0
      c(3) = 0.7876594617608471D0
      c(4) = 0.1D1


      end if
c
c

      aa(1) = 1.D0
      aa(2) = 1d0+c(2)
      aa(3) = 1d0+c(3)
      aa(4) = 1d0+c(4)


      vss(1,1,1) = 0.737280506621388D0
      vss(1,1,2) = 0.5702605427259742D0
      vss(1,1,3) = 0.4237909561944161D0
      vss(1,1,4) = 0.3494322101377552D0
      vss(1,2,1) = -0.8241873352699155D-1
      vss(1,2,2) = -0.4945467860500921D-1
      vss(1,2,3) = 0.4050873252770163D-1
      vss(1,2,4) = 0.113145853714089D0
      vss(1,3,1) = -0.1886811845882098D-1
      vss(1,3,2) = -0.1323354839475141D-1
      vss(1,3,3) = -0.4002796087543365D-3
      vss(1,3,4) = 0.8170306973806089D-2
      vss(1,4,1) = 0.3640063453644245D0
      vss(1,4,2) = 0.4924276842737864D0
      vss(1,4,3) = 0.5361005908866366D0
      vss(1,4,4) = 0.5292516291743497D0
      vss(2,1,1) = 0.17437746D-91
      vss(2,1,2) = 0.1379092013122435D-1
      vss(2,1,3) = 0.1758328714593773D-2
      vss(2,1,4) = 0.4505607022276855D-3
      vss(2,2,1) = 0.1D1
      vss(2,2,2) = 0.5162981615057852D0
      vss(2,2,3) = 0.3170279086894154D-1
      vss(2,2,4) = 0.6562999794283273D-2
      vss(2,3,1) = 0.4806015D-91
      vss(2,3,2) = 0.278273195396675D0
      vss(2,3,3) = 0.9119589977307768D0
      vss(2,3,4) = 0.967241218766513D0
      vss(2,4,1) = -0.13644872D-90
      vss(2,4,2) = 0.1916377229663154D0
      vss(2,4,3) = 0.5457988268568784D-1
      vss(2,4,4) = 0.2574522073697601D-1
      vss(3,1,1) = 0.1971532D-91
      vss(3,1,2) = 0.6831879755583199D-1
      vss(3,1,3) = 0.8767571397129687D1
      vss(3,1,4) = 0.425169763809326D2
      vss(3,2,1) = 0.8033644D-91
      vss(3,2,2) = 0.6804269103504088D0
      vss(3,2,3) = 0.7271550651590837D2
      vss(3,2,4) = 0.3302388932320554D3
      vss(3,3,1) = 0.1D1
      vss(3,3,2) = 0.2778510789827007D1
      vss(3,3,3) = 0.1635250440404749D3
      vss(3,3,4) = 0.628809673380971D3
      vss(3,4,1) = -0.15454554D-90
      vss(3,4,2) = -0.2527256497733248D1
      vss(3,4,3) = -0.2440081219535129D3
      vss(3,4,4) = -0.1000565542993959D4
      vss(4,1,1) = 0.2889016D-91
      vss(4,1,2) = 0.8900822227230133D1
      vss(4,1,3) = 0.2021745038925905D3
      vss(4,1,4) = 0.7094438857071275D3
      vss(4,2,1) = 0.11807857D-90
      vss(4,2,2) = 0.7377414176482106D2
      vss(4,2,3) = 0.1475652696154459D4
      vss(4,2,4) = 0.4941550595204215D4
      vss(4,3,1) = 0.80345D-91
      vss(4,3,2) = 0.1656249281630393D3
      vss(4,3,3) = 0.2439069600763088D4
      vss(4,3,4) = 0.7408577670496557D4
      vss(4,4,1) = 0.1D1
      vss(4,4,2) = -0.2472998921550905D3
      vss(4,4,3) = -0.4115896800810138D4
      vss(4,4,4) = -0.130585721514079D5
      wss(1,1,1) = 0.3559715077125859D-1
      wss(1,1,2) = 0.3271263703988755D-1
      wss(1,1,3) = 0.269277232182477D-1
      wss(1,1,4) = 0.2319084917324257D-1
      wss(1,2,1) = -0.1222393882782229D0
      wss(1,2,2) = -0.1785109611194919D0
      wss(1,2,3) = -0.2112378632597324D0
      wss(1,2,4) = -0.2200145802317251D0
      wss(1,3,1) = -0.1172773079710909D0
      wss(1,3,2) = -0.1602174452194821D0
      wss(1,3,3) = -0.1761777006239687D0
      wss(1,3,4) = -0.1749818625796969D0
      wss(1,4,1) = -0.2288944795547676D-1
      wss(1,4,2) = -0.3087653966701295D-1
      wss(1,4,3) = -0.3351934586961218D-1
      wss(1,4,4) = -0.3303501404894956D-1
      wss(2,1,1) = 0.1928104302096946D-92
      wss(2,1,2) = 0.1236157536238015D-2
      wss(2,1,3) = 0.162963446201919D-3
      wss(2,1,4) = 0.4232735863643932D-4
      wss(2,2,1) = 0.2063877941999665D-91
      wss(2,2,2) = 0.5292639657941568D-1
      wss(2,2,3) = 0.4600434968504916D-2
      wss(2,2,4) = 0.1052964093516254D-2
      wss(2,3,1) = 0.3794736889331792D-91
      wss(2,3,2) = -0.8810086609779139D-1
      wss(2,3,3) = -0.5570593048675531D-1
      wss(2,3,4) = 0.2881565099618327D-1
      wss(2,4,1) = 0.9027624995479898D-92
      wss(2,4,2) = -0.1116093808542199D-1
      wss(2,4,3) = -0.2944411790324903D-2
      wss(2,4,4) = -0.1266463017470281D-2
      wss(3,1,1) = 0.2179202900572329D-92
      wss(3,1,2) = 0.6617266324091891D-2
      wss(3,1,3) = 0.8651842594297778D0
      wss(3,1,4) = 0.4227028450091146D1
      wss(3,2,1) = 0.233528106436863D-91
      wss(3,2,2) = 0.1304049787542798D0
      wss(3,2,3) = 0.1517409229959125D2
      wss(3,2,4) = 0.7104321633643338D2
      wss(3,3,1) = 0.4297020482492912D-91
      wss(3,3,2) = 0.5921890488470773D0
      wss(3,3,3) = 0.5105056552602749D2
      wss(3,3,4) = 0.2197914614992998D3
      wss(3,4,1) = 0.1022588296724085D-91
      wss(3,4,2) = 0.4410925289488116D0
      wss(3,4,3) = 0.1974970813529859D2
      wss(3,4,4) = 0.7656938261827839D2
      wss(4,1,1) = 0.3191128889456822D-92
      wss(4,1,2) = 0.8783946011331707D0
      wss(4,1,3) = 0.2025196144741612D2
      wss(4,1,4) = 0.7149080353290433D2
      wss(4,2,1) = 0.3427510536692009D-91
      wss(4,2,2) = 0.1539951085997659D2
      wss(4,2,3) = 0.3267415532270721D3
      wss(4,2,4) = 0.1118149441254166D4
      wss(4,3,1) = 0.6317007028408447D-91
      wss(4,3,2) = 0.5176477020952794D2
      wss(4,3,3) = 0.9412338428105829D3
      wss(4,3,4) = 0.3066205067558872D4
      wss(4,4,1) = 0.1504391081489481D-91
      wss(4,4,2) = 0.200025750725675D2
      wss(4,4,3) = 0.3040699713546467D3
      wss(4,4,4) = 0.9458893110116933D3


      return

 500  continue

      an(1,1) = 0.2243471120859683D-2
      an(1,2) = -0.1036871827777886D-2
      an(1,3) = 0.7313354847739959D-3
      an(1,4) = -0.5056098356338289D-3
      an(1,5) = 0.1981196647206847D-3
      an(2,1) = 0.3108273511542618D-1
      an(2,2) = 0.8341706643245226D-2
      an(2,3) = -0.1571029336956402D-2
      an(2,4) = 0.706725930383004D-3
      an(2,5) = -0.2391112519788006D-3
      an(3,1) = 0.7589152043040584D-1
      an(3,2) = 0.8508291363322848D-1
      an(3,3) = 0.101255266506824D-1
      an(3,4) = -0.1000952110713696D-2
      an(3,5) = 0.1898877726666178D-3
      an(4,1) = 0.115333944047102D0
      an(4,2) = 0.1645006792668138D0
      an(4,3) = 0.8525394016214382D-1
      an(4,4) = 0.5009460194266192D-2
      an(4,5) = -0.9147817341039539D-4
      an(5,1) = 0.1355069134314881D0
      an(5,2) = 0.2034645680102714D0
      an(5,3) = 0.1298475476082324D0
      an(5,4) = 0.3118097095000808D-1
      an(5,5) = -0.12D-98
      bn(1) = 0.1355069134314881D0
      bn(2) = 0.2034645680102714D0
      bn(3) = 0.1298475476082324D0
      bn(4) = 0.3118097095000808D-1
      bn(5) = -0.7D-99
      c(1) = 0.5710419611451768D-1
      c(2) = 0.2768430136381238D0
      c(3) = 0.5835904323689168D0
      c(4) = 0.8602401356562194D0
      c(5) = 0.1D1
      b(1) = 0.1437135607912259D0
      b(2) = 0.2813560151494621D0
      b(3) = 0.3118265229757413D0
      b(4) = 0.2231039010835707D0
      b(5) = 0.4D-1

      bne(1) = 0.D0
      bne(2) = 0.6209328198572522D0
      bne(3) = -0.5082530842730005D0
      bne(4) = 0.686483734883042D0
      bne(5) = -0.2991634704672936D0


      vss(1,1,1) = 0.7371626170591345D0
      vss(1,1,2) = 0.6258043965510494D0
      vss(1,1,3) = 0.5050223878896125D0
      vss(1,1,4) = 0.4040848192716636D0
      vss(1,1,5) = 0.3570839619579335D0
      vss(1,2,1) = -0.8329661678083029D-1
      vss(1,2,2) = -0.6876182969625309D-1
      vss(1,2,3) = -0.1772057312364236D-1
      vss(1,2,4) = 0.5767966710439253D-1
      vss(1,2,5) = 0.1048165087465517D0
      vss(1,3,1) = -0.2077471237977784D-1
      vss(1,3,2) = -0.1797151994382234D-1
      vss(1,3,3) = -0.9289984732934516D-2
      vss(1,3,4) = 0.160803737968213D-2
      vss(1,3,5) = 0.7519529303094624D-2
      vss(1,4,1) = -0.6840786807026552D-2
      vss(1,4,2) = -0.5974260582663387D-2
      vss(1,4,3) = -0.3361040078129753D-2
      vss(1,4,4) = -0.18375604080198D-3
      vss(1,4,5) = 0.1495209519095813D-2
      vss(1,5,1) = 0.3737494989085002D0
      vss(1,5,2) = 0.4669032136716894D0
      vss(1,5,3) = 0.5253492100450941D0
      vss(1,5,4) = 0.5368112322850637D0
      vss(1,5,5) = 0.5290847904733243D0
      vss(2,1,1) = 0.171433734D-89
      vss(2,1,2) = 0.1087772423966831D-1
      vss(2,1,3) = 0.1209424390074593D-1
      vss(2,1,4) = 0.3375078135133181D-2
      vss(2,1,5) = 0.6172015599968827D-3
      vss(2,2,1) = 0.1D1
      vss(2,2,2) = 0.7860282057172373D0
      vss(2,2,3) = 0.334233634603881D0
      vss(2,2,4) = 0.628987454715301D-1
      vss(2,2,5) = 0.1007457615517635D-1
      vss(2,3,1) = 0.75731621D-89
      vss(2,3,2) = 0.9636059018470435D-1
      vss(2,3,3) = 0.4616118260603029D0
      vss(2,3,4) = 0.8498505906717781D0
      vss(2,3,5) = 0.96972069167906D0
      vss(2,4,1) = 0.45365885D-89
      vss(2,4,2) = 0.1411346012391125D-1
      vss(2,4,3) = 0.3896571189444054D-1
      vss(2,4,4) = 0.2394865769763749D-1
      vss(2,4,5) = 0.659146547707943D-2
      vss(2,5,1) = -0.2200336185D-88
      vss(2,5,2) = 0.9262001973447878D-1
      vss(2,5,3) = 0.1530945835406296D0
      vss(2,5,4) = 0.5992692802392111D-1
      vss(2,5,5) = 0.1299606512868731D-1
      vss(3,1,1) = 0.370721957D-89
      vss(3,1,2) = 0.3099970445533361D-2
      vss(3,1,3) = 0.8775924876318687D-3
      vss(3,1,4) = 0.5484669712192866D-1
      vss(3,1,5) = 0.6801247555733437D0
      vss(3,2,1) = 0.177326368D-88
      vss(3,2,2) = 0.3468393232635417D-1
      vss(3,2,3) = 0.8259354529318982D-2
      vss(3,2,4) = 0.4687371828898504D0
      vss(3,2,5) = 0.5603864517043714D1
      vss(3,3,1) = 0.1D1
      vss(3,3,2) = 0.3292953111347906D0
      vss(3,3,3) = 0.3382177689393231D-1
      vss(3,3,4) = 0.1381571390093793D1
      vss(3,3,5) = 0.1482310985750718D2
      vss(3,4,1) = 0.98924172D-89
      vss(3,4,2) = 0.4645025587872745D0
      vss(3,4,3) = 0.6878145544519901D0
      vss(3,4,4) = 0.4914725267195336D1
      vss(3,4,5) = 0.3842776800452244D2
      vss(3,5,1) = -0.478058218D-88
      vss(3,5,2) = 0.1684182273060474D0
      vss(3,5,3) = 0.2692267216371267D0
      vss(3,5,4) = -0.5819880537300909D1
      vss(3,5,5) = -0.5853486713464668D2
      vss(4,1,1) = -0.589942852D-89
      vss(4,1,2) = 0.8872315393018823D-1
      vss(4,1,3) = 0.2617924216994222D2
      vss(4,1,4) = 0.4616519039669968D3
      vss(4,1,5) = 0.1433851584085958D4
      vss(4,2,1) = -0.282237283D-88
      vss(4,2,2) = 0.7535988885088924D0
      vss(4,2,3) = 0.2017667251699052D3
      vss(4,2,4) = 0.3362761148867354D4
      vss(4,2,5) = 0.1021884862653052D5
      vss(4,3,1) = -0.262261347D-88
      vss(4,3,2) = 0.217983794588552D1
      vss(4,3,3) = 0.4461486932097213D3
      vss(4,3,4) = 0.6492123232249838D4
      vss(4,3,5) = 0.187818185746908D5
      vss(4,4,1) = 0.1D1
      vss(4,4,2) = 0.7303117089005176D1
      vss(4,4,3) = 0.7700745584424115D3
      vss(4,4,4) = 0.8794764649953346D4
      vss(4,4,5) = 0.2354853123678496D5
      vss(4,5,1) = 0.761009889D-88
      vss(4,5,2) = -0.9325277077329777D1
      vss(4,5,3) = -0.144316921899198D4
      vss(4,5,4) = -0.1911030093503753D5
      vss(4,5,5) = -0.5398205002209224D5
      vss(5,1,1) = -0.64320099D-89
      vss(5,1,2) = 0.6659648981263701D1
      vss(5,1,3) = 0.3893660743511819D3
      vss(5,1,4) = 0.4517819437597721D4
      vss(5,1,5) = 0.1240019450836504D5
      vss(5,2,1) = -0.307164787D-88
      vss(5,2,2) = 0.526941198298184D2
      vss(5,2,3) = 0.28456600131654D4
      vss(5,2,4) = 0.315150815775313D5
      vss(5,2,5) = 0.8494685868657433D5
      vss(5,3,1) = -0.284760683D-88
      vss(5,3,2) = 0.12470724492325D3
      vss(5,3,3) = 0.5535942493292079D4
      vss(5,3,4) = 0.5528530794608408D5
      vss(5,3,5) = 0.1434313133042538D6
      vss(5,4,1) = -0.170742491D-88
      vss(5,4,2) = 0.2477115568882784D3
      vss(5,4,3) = 0.7593580441316154D4
      vss(5,4,4) = 0.6469158651910924D5
      vss(5,4,5) = 0.1590543621874043D6
      vss(5,5,1) = 0.1D1
      vss(5,5,2) = -0.4307725706226105D3
      vss(5,5,3) = -0.1636354902212482D5
      vss(5,5,4) = -0.1560087954803223D6
      vss(5,5,5) = -0.3998317286865975D6
      wss(1,1,1) = 0.2269019844705136D-1
      wss(1,1,2) = 0.2181167215941908D-1
      wss(1,1,3) = 0.1939810646605529D-1
      wss(1,1,4) = 0.1656760200537627D-1
      wss(1,1,5) = 0.1504637326769702D-1
      wss(1,2,1) = -0.8873596803227226D-1
      wss(1,2,2) = -0.1173582697373427D0
      wss(1,2,3) = -0.1416098182258314D0
      wss(1,2,4) = -0.1552309697148278D0
      wss(1,2,5) = -0.1590072785149027D0
      wss(1,3,1) = -0.1117795052433773D0
      wss(1,3,2) = -0.1414460391124787D0
      wss(1,3,3) = -0.1615565585660109D0
      wss(1,3,4) = -0.1674522153880922D0
      wss(1,3,5) = -0.1662823944880138D0
      wss(1,4,1) = -0.827650145251378D-1
      wss(1,4,2) = -0.1036260298843751D0
      wss(1,4,3) = -0.1169003879139217D0
      wss(1,4,4) = -0.1197423044028455D0
      wss(1,4,5) = -0.1181685978051458D0
      wss(1,5,1) = -0.1498628427385529D-1
      wss(1,5,2) = -0.187079152538688D-1
      wss(1,5,3) = -0.2103222854362286D-1
      wss(1,5,4) = -0.2147440103309543D-1
      wss(1,5,5) = -0.2115678346109105D-1
      wss(2,1,1) = 0.1174849558482586D-90
      wss(2,1,2) = 0.6109269566219233D-3
      wss(2,1,3) = 0.7008192113784518D-3
      wss(2,1,4) = 0.1994954297722231D-3
      wss(2,1,5) = 0.3677675593957814D-4
      wss(2,2,1) = 0.159475077549197D-89
      wss(2,2,2) = 0.4240268747636307D-1
      wss(2,2,3) = 0.2861810433489201D-1
      wss(2,2,4) = 0.6492058403233217D-2
      wss(2,2,5) = 0.1106887725694985D-2
      wss(2,3,1) = 0.4326346022651087D-89
      wss(2,3,2) = -0.4262618900327403D-1
      wss(2,3,3) = -0.9611808667860523D-1
      wss(2,3,4) = -0.6802678733175451D-1
      wss(2,3,5) = -0.3092749199427303D-1
      wss(2,4,1) = 0.4468268400168761D-89
      wss(2,4,2) = -0.2186209617345075D-1
      wss(2,4,3) = -0.3736544863970047D-1
      wss(2,4,4) = -0.152664075982259D-1
      wss(2,4,5) = -0.3407684513452404D-2
      wss(2,5,1) = 0.9099562124009451D-90
      wss(2,5,2) = -0.3642781768491545D-2
      wss(2,5,3) = -0.5966385209197946D-2
      wss(2,5,4) = -0.2310150137404122D-2
      wss(2,5,5) = -0.4975627405184067D-3
      wss(3,1,1) = 0.2539322813347988D-90
      wss(3,1,2) = 0.1899688449570368D-3
      wss(3,1,3) = 0.5465020012346937D-4
      wss(3,1,4) = 0.3450033692090836D-2
      wss(3,1,5) = 0.4295596395375752D-1
      wss(3,2,1) = 0.3453290696471518D-89
      wss(3,2,2) = 0.4548774937627402D-2
      wss(3,2,3) = 0.11741414648231D-2
      wss(3,2,4) = 0.6969037023598544D-1
      wss(3,2,5) = 0.8474360781778998D0
      wss(3,3,1) = 0.9384235857469968D-89
      wss(3,3,2) = 0.3979506397615141D-1
      wss(3,3,3) = 0.673549241337138D-2
      wss(3,3,4) = 0.338856276274672D0
      wss(3,3,5) = 0.3902081659232218D1
      wss(3,4,1) = 0.9704216795471334D-89
      wss(3,4,2) = -0.6310041146459654D-1
      wss(3,4,3) = 0.4053037327603024D-1
      wss(3,4,4) = 0.8525465975500243D0
      wss(3,4,5) = 0.8379857774142792D1
      wss(3,5,1) = 0.1977330292254958D-89
      wss(3,5,2) = -0.605759873996739D-2
      wss(3,5,3) = -0.7263449497738283D-2
      wss(3,5,4) = 0.4770732148517634D0
      wss(3,5,5) = 0.327059233098903D1
      wss(4,1,1) = -0.4040771354402235D-90
      wss(4,1,2) = 0.5584734881616694D-2
      wss(4,1,3) = 0.1666366577452291D1
      wss(4,1,4) = 0.2959078471351698D2
      wss(4,1,5) = 0.9216443265451872D2
      wss(4,2,1) = -0.5495875798049633D-89
      wss(4,2,2) = 0.1123638648326945D0
      wss(4,2,3) = 0.3147265933539091D2
      wss(4,2,4) = 0.5384568053679358D3
      wss(4,2,5) = 0.165293253207162D4
      wss(4,3,1) = -0.1493675267259714D-88
      wss(4,3,2) = 0.5412066525028313D0
      wss(4,3,3) = 0.1324092467549608D3
      wss(4,3,4) = 0.2115490854448747D4
      wss(4,3,5) = 0.6334867778868976D4
      wss(4,4,1) = -0.1544747606907646D-88
      wss(4,4,2) = 0.1321405278532488D1
      wss(4,4,3) = 0.2319715861697163D3
      wss(4,4,4) = 0.3282731262878098D4
      wss(4,4,5) = 0.9456335879841202D4
      wss(4,5,1) = -0.3147702083628104D-89
      wss(4,5,2) = 0.6799125617821904D0
      wss(4,5,3) = 0.6833676856390137D2
      wss(4,5,4) = 0.8603240106168866D3
      wss(4,5,5) = 0.2400171700931751D4
      wss(5,1,1) = -0.4407116259566928D-90
      wss(5,1,2) = 0.4225803492514275D0
      wss(5,1,3) = 0.2494692146733649D2
      wss(5,1,4) = 0.2912116513361554D3
      wss(5,1,5) = 0.801237455509959D3
      wss(5,2,1) = -0.5986301032802897D-89
      wss(5,2,2) = 0.8119830443257868D1
      wss(5,2,3) = 0.4549548890528322D3
      wss(5,2,4) = 0.5148587101216173D4
      wss(5,2,5) = 0.1399479534264262D5
      wss(5,3,1) = -0.1625029037093647D-88
      wss(5,3,2) = 0.3535071446792359D2
      wss(5,3,3) = 0.1794331650659089D4
      wss(5,3,4) = 0.192729811848109D5
      wss(5,3,5) = 0.5138605481238074D5
      wss(5,4,1) = -0.1679130707501251D-88
      wss(5,4,2) = 0.6644374094265075D2
      wss(5,4,3) = 0.2801823057678908D4
      wss(5,4,4) = 0.2779119984240949D5
      wss(5,4,5) = 0.7212858609936881D5
      wss(5,5,1) = -0.3420237737556848D-89
      wss(5,5,2) = 0.2125517021298784D2
      wss(5,5,3) = 0.7382492576807989D3
      wss(5,5,4) = 0.686766985079168D4
      wss(5,5,5) = 0.1747602236753508D5

      aa(1) = 1D0
      aa(2) = 1d0+c(2)
      aa(3) = 1d0+c(3)
      aa(4) = 1d0+c(4)
      aa(5) = 1d0+c(5)


      return


 600  continue

      c(1) = 0.3980985705146874D-1
      c(2) = 0.1980134178736082D0
      c(3) = 0.4379748102473861D0
      c(4) = 0.6954642733536361D0
      c(5) = 0.9014649142011736D0
      c(6) = 0.1D1
      an(1,1) = 0.1274167744433679D-2
      an(1,2) = -0.1095400226586424D-2
      an(1,3) = 0.1344465624980305D-2
      an(1,4) = -0.1464644255096784D-2
      an(1,5) = 0.1292565115958251D-2
      an(1,6) = -0.5587416444598391D-3
      an(2,1) = 0.1547059341883339D-1
      an(2,2) = 0.5293470140520667D-2
      an(2,3) = -0.2191834445963264D-2
      an(2,4) = 0.2000523479449754D-2
      an(2,5) = -0.1691846489382298D-2
      an(2,6) = 0.7237507255358312D-3
      an(3,1) = 0.4047129670403428D-1
      an(3,2) = 0.4861718370803604D-1
      an(3,3) = 0.8205095861205075D-2
      an(3,4) = -0.2346719041383018D-2
      an(3,5) = 0.1638922671984232D-2
      an(3,6) = -0.674812698259661D-3
      an(4,1) = 0.658684618079735D-1
      an(4,2) = 0.1044617594252577D0
      an(4,3) = 0.6536592191132361D-1
      an(4,4) = 0.7059790594744501D-2
      an(4,5) = -0.1397990400879028D-2
      an(4,6) = 0.4773344172302175D-3
      an(5,1) = 0.8694484453428965D-1
      an(5,2) = 0.1463139747658316D0
      an(5,3) = 0.1213277037080849D0
      an(5,4) = 0.4898879199164998D-1
      an(5,5) = 0.2960864088848145D-2
      an(5,6) = -0.2166833208397298D-3
      an(6,1) = 0.9678159022665168D-1
      an(6,2) = 0.1671746380943696D0
      an(6,3) = 0.1463869870846698D0
      an(6,4) = 0.7390887007261667D-1
      an(6,5) = 0.1574791452169228D-1
      an(6,6) = 0.35805D-98
      bn(1) = 0.9678159022665168D-1
      bn(2) = 0.1671746380943696D0
      bn(3) = 0.1463869870846698D0
      bn(4) = 0.7390887007261667D-1
      bn(5) = 0.1574791452169228D-1
      bn(6) = 0.35805D-98
      b(1) = 0.1007941926267404D0
      b(2) = 0.2084506671559539D0
      b(3) = 0.2604633915947875D0
      b(4) = 0.242693594234485D0
      b(5) = 0.1598203766102555D0
      b(6) = 0.2777777777777778D-1
      aa(1) = 1.D0
      aa(2) = 0.1198013417873608D1
      aa(3) = 0.1437974810247386D1
      aa(4) = 0.1695464273353636D1
      aa(5) = 0.1901464914201174D1
      aa(6) = 0.2D1
      be(1) = 0.D0
      be(2) = 0.5317251316291001D0
      be(3) = -0.2769646937792649D0
      be(4) = 0.8964094447073625D0
      be(5) = -0.4441462412265267D0
      be(6) = 0.292976358669329D0
      bne(1) = 0.D0
      bne(2) = 0.4878082266289541D0
      bne(3) = -0.4124431687270232D0
      bne(4) = 0.7873227228597435D0
      bne(5) = -0.6682570567168819D0
      bne(6) = 0.3055692759552075D0

      vss(1,1,1) = 0.7371210656220225D0
      vss(1,1,2) = 0.6581335938879463D0
      vss(1,1,3) = 0.5630010814224243D0
      vss(1,1,4) = 0.465680952894569D0
      vss(1,1,5) = 0.3934737526617886D0
      vss(1,1,6) = 0.3610651733187308D0
      vss(1,2,1) = -0.8359104894408238D-1
      vss(1,2,2) = -0.7631449360904815D-1
      vss(1,2,3) = -0.4733768793121524D-1
      vss(1,2,4) = 0.7667218746322173D-2
      vss(1,2,5) = 0.6757317838444819D-1
      vss(1,2,6) = 0.1005196027072396D0
      vss(1,3,1) = -0.2134981670891672D-1
      vss(1,3,2) = -0.1988155999566891D-1
      vss(1,3,3) = -0.1454026464940092D-1
      vss(1,3,4) = -0.560180627556167D-2
      vss(1,3,5) = 0.2867794116636934D-2
      vss(1,3,6) = 0.7052806351790091D-2
      vss(1,4,1) = -0.7889589182873309D-2
      vss(1,4,2) = -0.7379333341294802D-2
      vss(1,4,3) = -0.5560898657842525D-2
      vss(1,4,4) = -0.2598593646361707D-2
      vss(1,4,5) = 0.1313702308904538D-3
      vss(1,4,6) = 0.1453178992028402D-2
      vss(1,5,1) = -0.3092290021338287D-2
      vss(1,5,2) = -0.2896790610965097D-2
      vss(1,5,3) = -0.2205135980023655D-2
      vss(1,5,4) = -0.1088872493689964D-2
      vss(1,5,5) = -0.6980175755240466D-4
      vss(1,5,6) = 0.4203059734829264D-3
      vss(1,6,1) = 0.3788016792351882D0
      vss(1,6,2) = 0.4483385836690307D0
      vss(1,6,3) = 0.506642905796058D0
      vss(1,6,4) = 0.5359411007747222D0
      vss(1,6,5) = 0.5360237063637882D0
      vss(1,6,6) = 0.5294889326567281D0
      vss(2,1,1) = 0.8630761688D-88
      vss(2,1,2) = 0.696209982559561D-2
      vss(2,1,3) = 0.1411245666030413D-1
      vss(2,1,4) = 0.1018165725356068D-1
      vss(2,1,5) = 0.4078792197653709D-2
      vss(2,1,6) = 0.187562150556413D-2
      vss(2,2,1) = 0.1D1
      vss(2,2,2) = 0.8921931105205456D0
      vss(2,2,3) = 0.5828399034285181D0
      vss(2,2,4) = 0.2477536744449984D0
      vss(2,2,5) = 0.772878924176651D-1
      vss(2,2,6) = 0.3246295936026519D-1
      vss(2,3,1) = 0.4953935494D-87
      vss(2,3,2) = 0.4223059466601683D-1
      vss(2,3,3) = 0.2336982125704588D0
      vss(2,3,4) = 0.5659322820400471D0
      vss(2,3,5) = 0.8237768717938019D0
      vss(2,3,6) = 0.9151817790834255D0
      vss(2,4,1) = 0.3919393342D-87
      vss(2,4,2) = 0.7482737119191608D-2
      vss(2,4,3) = 0.3012823449510456D-1
      vss(2,4,4) = 0.4243730507208818D-1
      vss(2,4,5) = 0.2842462421500857D-1
      vss(2,4,6) = 0.1669828380414069D-1
      vss(2,5,1) = 0.2273222386D-87
      vss(2,5,2) = 0.233458646340077D-2
      vss(2,5,3) = 0.8771688703226713D-2
      vss(2,5,4) = 0.1127367035948414D-1
      vss(2,5,5) = 0.6897482421725558D-2
      vss(2,5,6) = 0.385318717406102D-2
      vss(2,6,1) = -0.16529726512D-86
      vss(2,6,2) = 0.4879687140524962D-1
      vss(2,6,3) = 0.1304495041423877D0
      vss(2,6,4) = 0.1224214108298215D0
      vss(2,6,5) = 0.5953433695414519D-1
      vss(2,6,6) = 0.2992816907254349D-1
      vss(3,1,1) = 0.60951465596D-86
      vss(3,1,2) = 0.3153134888889783D-2
      vss(3,1,3) = 0.1488406222269677D-2
      vss(3,1,4) = 0.5929110463824275D-3
      vss(3,1,5) = 0.1010266669325444D-2
      vss(3,1,6) = 0.1952495911028695D-3
      vss(3,2,1) = 0.317940152661D-85
      vss(3,2,2) = 0.3789684193890523D-1
      vss(3,2,3) = 0.1531520787838038D-1
      vss(3,2,4) = 0.5484154383915141D-2
      vss(3,2,5) = 0.8796990676262922D-2
      vss(3,2,6) = 0.166038116372924D-2
      vss(3,3,1) = 0.1D1
      vss(3,3,2) = 0.6843022449639048D0
      vss(3,3,3) = 0.9902101459174182D-1
      vss(3,3,4) = 0.2241237189315692D-1
      vss(3,3,5) = 0.2910950763697966D-1
      vss(3,3,6) = 0.5094960429272997D-2
      vss(3,4,1) = 0.27264711473D-85
      vss(3,4,2) = 0.1644333115567796D0
      vss(3,4,3) = 0.7879452003626522D0
      vss(3,4,4) = 0.8744928674771684D0
      vss(3,4,5) = 0.2246273725265967D0
      vss(3,4,6) = 0.2801470771465231D-1
      vss(3,5,1) = 0.15757790498D-85
      vss(3,5,2) = 0.2222708290932128D-1
      vss(3,5,3) = 0.3030621458072422D-1
      vss(3,5,4) = 0.5060996112288214D-1
      vss(3,5,5) = 0.5916314710419657D0
      vss(3,5,6) = 0.9232263836767512D0
      vss(3,6,1) = -0.1155543213473D-84
      vss(3,6,2) = 0.8798738374219928D-1
      vss(3,6,3) = 0.659239563642317D-1
      vss(3,6,4) = 0.4640773407649502D-1
      vss(3,6,5) = 0.1448243914488696D0
      vss(3,6,6) = 0.4180831742449142D-1
      vss(4,1,1) = 0.42297372447D-86
      vss(4,1,2) = 0.1003620925785823D-2
      vss(4,1,3) = 0.5099072567655288D-8
      vss(4,1,4) = 0.275995246584966D1
      vss(4,1,5) = 0.7210633268453077D2
      vss(4,1,6) = 0.2420565486189324D3
      vss(4,2,1) = 0.22053988651D-85
      vss(4,2,2) = 0.8737179992123907D-2
      vss(4,2,3) = 0.4104561136698795D-7
      vss(4,2,4) = 0.2100247020521935D2
      vss(4,2,5) = 0.5310278348868424D3
      vss(4,2,6) = 0.1759559300351788D4
      vss(4,3,1) = 0.2401479223D-85
      vss(4,3,2) = 0.2889060577288278D-1
      vss(4,3,3) = 0.10715097872892D-6
      vss(4,3,4) = 0.4724227178201099D2
      vss(4,3,5) = 0.1102464706038328D4
      vss(4,3,6) = 0.3542731321822363D4
      vss(4,4,1) = 0.1D1
      vss(4,4,2) = 0.2221330671005495D0
      vss(4,4,3) = 0.3408413375016119D-6
      vss(4,4,4) = 0.1036648737009846D3
      vss(4,4,5) = 0.2050791338021851D4
      vss(4,4,6) = 0.621806485255626D4
      vss(4,5,1) = 0.10913151974D-85
      vss(4,5,2) = 0.5949167035343624D0
      vss(4,5,3) = 0.2142874883833451D-5
      vss(4,5,4) = 0.2072303356794601D3
      vss(4,5,5) = 0.2997501288712887D4
      vss(4,5,6) = 0.8271873873879244D4
      vss(4,6,1) = -0.80101224493D-85
      vss(4,6,2) = 0.1443188226742956D0
      vss(4,6,3) = 0.999997362988116D0
      vss(4,6,4) = -0.3808999038335247D3
      vss(4,6,5) = -0.6752891500344439D4
      vss(4,6,6) = -0.2003328589722859D5
      vss(5,1,1) = -0.2099182781D-86
      vss(5,1,2) = 0.8004140898058004D-1
      vss(5,1,3) = 0.4336935328784321D2
      vss(5,1,4) = 0.1856949965998174D4
      vss(5,1,5) = 0.1694623749860255D5
      vss(5,1,6) = 0.4196739746584531D5
      vss(5,2,1) = -0.1094483213D-85
      vss(5,2,2) = 0.6266949016022258D0
      vss(5,2,3) = 0.3211133867583405D3
      vss(5,2,4) = 0.1320129552192196D5
      vss(5,2,5) = 0.117624837060652D6
      vss(5,2,6) = 0.2885263610379565D6
      vss(5,3,1) = -0.11917131607D-85
      vss(5,3,2) = 0.1517294148624002D1
      vss(5,3,3) = 0.6752719581184592D3
      vss(5,3,4) = 0.2526076412466838D5
      vss(5,3,5) = 0.2135262596971442D6
      vss(5,3,6) = 0.5131395409358237D6
      vss(5,4,1) = -0.9373050434D-86
      vss(5,4,2) = 0.395199624107312D1
      vss(5,4,3) = 0.1288076479854708D4
      vss(5,4,4) = 0.4047718942634008D5
      vss(5,4,5) = 0.3133815446245724D6
      vss(5,4,6) = 0.7289829288439548D6
      vss(5,5,1) = 0.1D1
      vss(5,5,2) = 0.1208014106207161D2
      vss(5,5,3) = 0.1965393733384294D4
      vss(5,5,4) = 0.4708222683112834D5
      vss(5,5,5) = 0.3245043799852836D6
      vss(5,5,6) = 0.725094046412918D6
      vss(5,6,1) = 0.39749004166D-85
      vss(5,6,2) = -0.1725616776235153D2
      vss(5,6,3) = -0.4292224911403645D4
      vss(5,6,4) = -0.1278774258700569D6
      vss(5,6,5) = -0.9859822588662547D6
      vss(5,6,6) = -0.2297709274696498D7
      vss(6,1,1) = -0.20173923686D-86
      vss(6,1,2) = 0.4606442304729594D1
      vss(6,1,3) = 0.5112530871815369D3
      vss(6,1,4) = 0.1385972175847366D5
      vss(6,1,5) = 0.105302541545596D6
      vss(6,1,6) = 0.2450543829896044D6
      vss(6,2,1) = -0.10496750199D-85
      vss(6,2,2) = 0.3488515212763787D2
      vss(6,2,3) = 0.3686242117960601D4
      vss(6,2,4) = 0.964076681049273D5
      vss(6,2,5) = 0.7171244928166298D6
      vss(6,2,6) = 0.1654761428564075D7
      vss(6,3,1) = -0.11394729417D-85
      vss(6,3,2) = 0.775268267110139D2
      vss(6,3,3) = 0.7283408691357278D4
      vss(6,3,4) = 0.1758243383131631D6
      vss(6,3,5) = 0.1250105203216348D7
      vss(6,3,6) = 0.2834102010213437D7
      vss(6,4,1) = -0.8937676294D-86
      vss(6,4,2) = 0.1657134498337557D3
      vss(6,4,3) = 0.1235034003881019D5
      vss(6,4,4) = 0.2599888179564531D6
      vss(6,4,5) = 0.172174245941719D7
      vss(6,4,6) = 0.3800536991943994D7
      vss(6,5,1) = -0.5153761895D-86
      vss(6,5,2) = 0.3135520681235324D3
      vss(6,5,3) = 0.155876019854145D5
      vss(6,5,4) = 0.2717770674781887D6
      vss(6,5,5) = 0.1650403944966554D7
      vss(6,5,6) = 0.3532380941274346D7
      vss(6,6,1) = 0.1D1
      vss(6,6,2) = -0.5952839391006694D3
      vss(6,6,3) = -0.3941784592072411D5
      vss(6,6,4) = -0.8178566136112058D6
      vss(6,6,5) = -0.5444677641962318D7
      vss(6,6,6) = -0.1206683475498546D8
      wss(1,1,1) = 0.1572337331305659D-1
      wss(1,1,2) = 0.154025989269822D-1
      wss(1,1,3) = 0.1435136756142564D-1
      wss(1,1,4) = 0.127315806919054D-1
      wss(1,1,5) = 0.1125119137955516D-1
      wss(1,1,6) = 0.1052027432361482D-1
      wss(1,2,1) = -0.6577790409883515D-1
      wss(1,2,2) = -0.8124329553780172D-1
      wss(1,2,3) = -0.9702678916138957D-1
      wss(1,2,4) = -0.1094595579012339D0
      wss(1,2,5) = -0.1157584755660238D0
      wss(1,2,6) = -0.1176194144354173D0
      wss(1,3,1) = -0.9352241142115749D-1
      wss(1,3,2) = -0.1119721695802647D0
      wss(1,3,3) = -0.1283693483391501D0
      wss(1,3,4) = -0.1379878412916789D0
      wss(1,3,5) = -0.139855049329321D0
      wss(1,3,6) = -0.1390517555231841D0
      wss(1,4,1) = -0.9036777455432457D-1
      wss(1,4,2) = -0.1073437305749294D0
      wss(1,4,3) = -0.1218478560414663D0
      wss(1,4,4) = -0.1295317554009355D0
      wss(1,4,5) = -0.1300771680825639D0
      wss(1,4,6) = -0.1287445351459449D0
      wss(1,5,1) = -0.6033535472846147D-1
      wss(1,5,2) = -0.7146168958045767D-1
      wss(1,5,3) = -0.8082558439954256D-1
      wss(1,5,4) = -0.8558167691213312D-1
      wss(1,5,5) = -0.8566204429256234D-1
      wss(1,5,6) = -0.8464993383935949D-1
      wss(1,6,1) = -0.1053440110013876D-1
      wss(1,6,2) = -0.1246522043231289D-1
      wss(1,6,3) = -0.1408209768280782D-1
      wss(1,6,4) = -0.1489161236041915D-1
      wss(1,6,5) = -0.1488997110947396D-1
      wss(1,6,6) = -0.1470656071256203D-1
      wss(2,1,1) = 0.4034799355185552D-89
      wss(2,1,2) = 0.2680139351751601D-3
      wss(2,1,3) = 0.5583964933949513D-3
      wss(2,1,4) = 0.4115494113902575D-3
      wss(2,1,5) = 0.167046416669831D-3
      wss(2,1,6) = 0.7722971450215912D-4
      wss(2,2,1) = 0.6244272372151547D-88
      wss(2,2,2) = 0.2695779579475722D-1
      wss(2,2,3) = 0.3047928460552199D-1
      wss(2,2,4) = 0.166799884117547D-1
      wss(2,2,5) = 0.5857827703574237D-2
      wss(2,2,6) = 0.2567877882809177D-2
      wss(2,3,1) = 0.2058003615323045D-87
      wss(2,3,2) = -0.188406555448267D-1
      wss(2,3,3) = -0.6172608470739497D-1
      wss(2,3,4) = -0.8038661516785319D-1
      wss(2,3,5) = -0.6088008807457981D-1
      wss(2,3,6) = -0.4337509436511458D-1
      wss(2,4,1) = 0.3013853328534143D-87
      wss(2,4,2) = -0.1313322710615662D-1
      wss(2,4,3) = -0.366428002273389D-1
      wss(2,4,4) = -0.3635542594791507D-1
      wss(2,4,5) = -0.1866600977632936D-1
      wss(2,4,6) = -0.9668412731700345D-2
      wss(2,5,1) = 0.248382088597037D-87
      wss(2,5,2) = -0.7949799828939233D-2
      wss(2,5,3) = -0.2141237610738032D-1
      wss(2,5,4) = -0.2028416554594784D-1
      wss(2,5,5) = -0.9951113226272468D-2
      wss(2,5,6) = -0.5025856165710851D-2
      wss(2,6,1) = 0.4695059607769272D-88
      wss(2,6,2) = -0.1346935145824845D-2
      wss(2,6,3) = -0.359214917747644D-2
      wss(2,6,4) = -0.3361149867253606D-2
      wss(2,6,5) = -0.1630141241977476D-2
      wss(2,6,6) = -0.8183180906248786D-3
      wss(3,1,1) = 0.2851483914903219D-87
      wss(3,1,2) = 0.1332289604217836D-3
      wss(3,1,3) = 0.6374603731055724D-4
      wss(3,1,4) = 0.2565841509958018D-4
      wss(3,1,5) = 0.4399739874378247D-4
      wss(3,1,6) = 0.8525139055986276D-5
      wss(3,2,1) = 0.4400361184627248D-86
      wss(3,2,2) = 0.3557447491995136D-2
      wss(3,2,3) = 0.1545779977681442D-2
      wss(3,2,4) = 0.5816798498750469D-3
      wss(3,2,5) = 0.9596236432008112D-3
      wss(3,2,6) = 0.1831267195535253D-3
      wss(3,3,1) = 0.1445977428268324D-85
      wss(3,3,2) = 0.4795214999161297D-1
      wss(3,3,3) = 0.124508132912456D-1
      wss(3,3,4) = 0.3721225836966168D-2
      wss(3,3,5) = 0.5520887165684987D-2
      wss(3,3,6) = 0.1014317340428267D-2
      wss(3,4,1) = 0.2111915727239902D-85
      wss(3,4,2) = -0.3874930429856453D-1
      wss(3,4,3) = -0.5787022430576256D-1
      wss(3,4,4) = 0.3828601785253039D-1
      wss(3,4,5) = 0.2525532528080399D-1
      wss(3,4,6) = 0.3916434660952839D-2
      wss(3,5,1) = 0.1737344066252576D-85
      wss(3,5,2) = -0.1541017685343399D-1
      wss(3,5,3) = -0.1227393489432845D-1
      wss(3,5,4) = -0.9959288938081319D-2
      wss(3,5,5) = -0.4431583252989078D-1
      wss(3,5,6) = -0.243079870513128D-1
      wss(3,6,1) = 0.3281458541948916D-86
      wss(3,6,2) = -0.2376687680841134D-2
      wss(3,6,3) = -0.1752443884994967D-2
      wss(3,6,4) = -0.1196025064718662D-2
      wss(3,6,5) = -0.3541136383706902D-2
      wss(3,6,6) = -0.972776908236135D-3
      wss(4,1,1) = 0.197894639342239D-87
      wss(4,1,2) = 0.4370902215296425D-4
      wss(4,1,3) = 0.2240254972210779D-9
      wss(4,1,4) = 0.1220717436464756D0
      wss(4,1,5) = 0.3202251422256092D1
      wss(4,1,6) = 0.1076757765775937D2
      wss(4,2,1) = 0.3052934172537964D-86
      wss(4,2,2) = 0.9531976311457186D-3
      wss(4,2,3) = 0.464391053670695D-8
      wss(4,2,4) = 0.2439022483959466D1
      wss(4,2,5) = 0.6261271013451971D2
      wss(4,2,6) = 0.2087244033869846D3
      wss(4,3,1) = 0.1002885535806031D-85
      wss(4,3,2) = 0.5481914433269192D-2
      wss(4,3,3) = 0.2371175566815254D-7
      wss(4,3,4) = 0.1155337885515963D2
      wss(4,3,5) = 0.2848389645594181D3
      wss(4,3,6) = 0.9349607915138729D3
      wss(4,4,1) = 0.1464337467958678D-85
      wss(4,4,2) = 0.2503170240305049D-1
      wss(4,4,3) = 0.6962285515842134D-7
      wss(4,4,4) = 0.2816963372669151D2
      wss(4,4,5) = 0.6393624837315335D3
      wss(4,4,6) = 0.2038450840715898D4
      wss(4,5,1) = 0.1204384786312191D-85
      wss(4,5,2) = -0.4429180468636767D-1
      wss(4,5,3) = 0.1887281932885638D-6
      wss(4,5,4) = 0.4305452046986715D2
      wss(4,5,5) = 0.8355481183202024D3
      wss(4,5,6) = 0.2541390323114083D4
      wss(4,6,1) = 0.2274624212404132D-86
      wss(4,6,2) = -0.3527594009646273D-2
      wss(4,6,3) = 0.6023276416721494D-4
      wss(4,6,4) = 0.1298462262064298D2
      wss(4,6,5) = 0.2122526163430612D3
      wss(4,6,6) = 0.619093433591818D3
      wss(5,1,1) = -0.9821391852552257D-88
      wss(5,1,2) = 0.3528051930031409D-2
      wss(5,1,3) = 0.1924736195121182D1
      wss(5,1,4) = 0.8284235970975967D2
      wss(5,1,5) = 0.7584047904490301D3
      wss(5,1,6) = 0.1880610347856015D4
      wss(5,2,1) = -0.1515118693976185D-86
      wss(5,2,2) = 0.7182235838051769D-1
      wss(5,2,3) = 0.3776782494732849D2
      wss(5,2,4) = 0.1582236904244343D4
      wss(5,2,5) = 0.1425519067696325D5
      wss(5,2,6) = 0.3512247962609432D5
      wss(5,3,1) = -0.4976998625021486D-86
      wss(5,3,2) = 0.3530677704136242D0
      wss(5,3,3) = 0.1729298727365131D3
      wss(5,3,4) = 0.6907737406252564D4
      wss(5,3,5) = 0.6060186425007903D5
      wss(5,3,6) = 0.147775245210332D6
      wss(5,4,1) = -0.7266796749840208D-86
      wss(5,4,2) = 0.937990022160153D0
      wss(5,4,3) = 0.3930786784040078D3
      wss(5,4,4) = 0.1438897521836706D5
      wss(5,4,5) = 0.1208008184853791D6
      wss(5,4,6) = 0.2898008480732303D6
      wss(5,5,1) = -0.5976625065193008D-86
      wss(5,5,2) = 0.1772845946966115D1
      wss(5,5,3) = 0.524858610514094D3
      wss(5,5,4) = 0.1677394623913992D5
      wss(5,5,5) = 0.1328646730825584D6
      wss(5,5,6) = 0.3123904760004168D6
      wss(5,6,1) = -0.1128743614072408D-86
      wss(5,6,2) = 0.7699545390404183D0
      wss(5,6,3) = 0.1360799391633608D3
      wss(5,6,4) = 0.3872189420527663D4
      wss(5,6,5) = 0.2941762132172783D5
      wss(5,6,6) = 0.6823870916776926D5
      wss(6,1,1) = -0.9442225715048233D-88
      wss(6,1,2) = 0.2038618529208165D0
      wss(6,1,3) = 0.2276626435447276D2
      wss(6,1,4) = 0.6200941717228755D3
      wss(6,1,5) = 0.4724821019256816D4
      wss(6,1,6) = 0.110081458348921D5
      wss(6,2,1) = -0.1454485740319482D-86
      wss(6,2,2) = 0.4060295026512563D1
      wss(6,2,3) = 0.438929561690765D3
      wss(6,2,4) = 0.1167222513371561D5
      wss(6,2,5) = 0.8768066752446977D5
      wss(6,2,6) = 0.2031199397008444D6
      wss(6,3,1) = -0.4770510543772397D-86
      wss(6,3,2) = 0.1911632381954917D2
      wss(6,3,3) = 0.1947535736711314D4
      wss(6,3,4) = 0.497375421714987D5
      wss(6,3,5) = 0.3651962174159741D6
      wss(6,3,6) = 0.8384912206453585D6
      wss(6,4,1) = -0.6955688523655979D-86
      wss(6,4,2) = 0.4600143139931427D2
      wss(6,4,3) = 0.4173425505451332D4
      wss(6,4,4) = 0.995172866390151D5
      wss(6,4,5) = 0.7051449069206952D6
      wss(6,4,6) = 0.1597498381529067D7
      wss(6,5,1) = -0.5715367097867036D-86
      wss(6,5,2) = 0.6840232760565543D2
      wss(6,5,3) = 0.5067972895245975D4
      wss(6,5,4) = 0.1099752880292815D6
      wss(6,5,5) = 0.7461783527610653D6
      wss(6,5,6) = 0.1664559047156363D7
      wss(6,6,1) = -0.1078964105593898D-86
      wss(6,6,2) = 0.1993121522878036D2
      wss(6,6,3) = 0.120788244449693D4
      wss(6,6,4) = 0.244288046737571D5
      wss(6,6,5) = 0.1610480976670908D6
      wss(6,6,6) = 0.3557745841308963D6

      return
      end
c-----!-----------------------------------------------------------------
c --- S U B R O U T I N E   R E C O E F F
c-----!-----------------------------------------------------------------
c
c Recompute the RK-coefficients for a stepsize change
c hnew=alpha*h.
c
      subroutine recoeff (s,aa,c,alpha,vss,wss,tmp,y,xinv)
      implicit none

      integer s
      real*8 aa(s),c(s)
      real*8 alpha,vss(s,s,s),wss(s,s,s),xinv(2*s,2*s)
      real*8 tmp(s,2*s),y(s,2*s)


      integer i,ii,j,i1,i2

      if (s.eq.3) goto 300
      if (s.eq.4) goto 400
      if (s.eq.5) goto 500


      stop "vs nur fuer s=3,..,5 implementiert!"

 300  continue
      xinv(1,1) = -0.1067102363878277D3
      xinv(1,2) = -0.2414330706801914D3
      xinv(1,3) = 0.3491433070680191D3
      xinv(1,4) = -0.2602041028867288D2
      xinv(1,5) = -0.1254816342116934D3
      xinv(1,6) = -0.4292925127559876D2
      xinv(2,1) = 0.3691837967430268D3
      xinv(2,2) = 0.9174456685847274D3
      xinv(2,3) = -0.1286629465327754D4
      xinv(2,4) = 0.8367755909695693D2
      xinv(2,5) = 0.4527278975326853D3
      xinv(2,6) = 0.1595183657883066D3
      xinv(3,1) = -0.4899675309990711D3
      xinv(3,2) = -0.1364096849343082D4
      xinv(3,3) = 0.1854064380342153D4
      xinv(3,4) = -0.1054153181310017D3
      xinv(3,5) = -0.6366642958808186D3
      xinv(3,6) = -0.2322119025302268D3
      xinv(4,1) = 0.3154878892152585D3
      xinv(4,2) = 0.9898755897887848D3
      xinv(4,3) = -0.1305363479004043D4
      xinv(4,4) = 0.6528368218355879D2
      xinv(4,5) = 0.4361421847347563D3
      xinv(4,6) = 0.1655183657883066D3
      xinv(5,1) = -0.9925719046627242D2
      xinv(5,2) = -0.3500779524862776D3
      xinv(5,3) = 0.44933514295255D3
      xinv(5,4) = -0.1992959491857567D2
      xinv(5,5) = -0.1457949008993308D3
      xinv(5,6) = -0.5782823076116511D2
      xinv(6,1) = 0.1226327189488592D2
      xinv(6,2) = 0.4828661413603829D2
      xinv(6,3) = -0.605498860309242D2
      xinv(6,4) = 0.2404082057734575D1
      xinv(6,5) = 0.1907074872440124D2
      xinv(6,6) = 0.7932652990377571D1

      goto 1000

 400  continue
      xinv(1,1) = -0.2056463062270675D4
      xinv(1,2) = -0.6638213150793579D4
      xinv(1,3) = -0.3663951489975653D4
      xinv(1,4) = 0.1235962770303991D5
      xinv(1,5) = -0.24413205201554D3
      xinv(1,6) = -0.2154429709669696D4
      xinv(1,7) = -0.3550468467886089D4
      xinv(1,8) = -0.8075230906118085D3
      xinv(2,1) = 0.9575862857597645D4
      xinv(2,2) = 0.3313556426656857D5
      xinv(2,3) = 0.1935149999034988D5
      xinv(2,4) = -0.620629271145161D5
      xinv(2,5) = 0.1107812064422838D4
      xinv(2,6) = 0.1040216711929215D5
      xinv(2,7) = 0.176755347086502D5
      xinv(2,8) = 0.4068105553780605D4
      xinv(3,1) = -0.1881522580326026D5
      xinv(3,2) = -0.6979664351781238D5
      xinv(3,3) = -0.4331591020176612D5
      xinv(3,4) = 0.1319277795228388D6
      xinv(3,5) = -0.2131109374720695D4
      xinv(3,6) = -0.2120279471315534D5
      xinv(3,7) = -0.3721684895084965D5
      xinv(3,8) = -0.867942719918242D4
      xinv(4,1) = 0.2026301074207982D5
      xinv(4,2) = 0.8041746831558644D5
      xinv(4,3) = 0.5324312367360044D5
      xinv(4,4) = -0.1539236027312667D6
      xinv(4,5) = 0.2255051572066687D4
      xinv(4,6) = 0.2366409884437793D5
      xinv(4,7) = 0.4297018783166164D5
      xinv(4,8) = 0.1016852575050078D5
      xinv(5,1) = -0.129379323560964D5
      xinv(5,2) = -0.5474923041161093D5
      xinv(5,3) = -0.3879799668097592D5
      xinv(5,4) = 0.1064851594486833D6
      xinv(5,5) = -0.1418675487916385D4
      xinv(5,6) = -0.1563070777839437D5
      xinv(5,7) = -0.293907342776322D5
      xinv(5,8) = -0.706743701441803D4
      xinv(6,1) = 0.4903794085392577D4
      xinv(6,2) = 0.2203736600407478D5
      xinv(6,3) = 0.1675421078707778D5
      xinv(6,4) = -0.4369537087654513D5
      xinv(6,5) = 0.5309793360435537D3
      xinv(6,6) = 0.6115617601287364D4
      xinv(6,7) = 0.1191337546346276D5
      xinv(6,8) = 0.2915214666752651D4
      xinv(7,1) = -0.1022640764205055D4
      xinv(7,2) = -0.4859589990456991D4
      xinv(7,3) = -0.3968714729068658D4
      xinv(7,4) = 0.9850945483730704D4
      xinv(7,5) = -0.109539645429474D3
      xinv(7,6) = -0.1313528427191663D4
      xinv(7,7) = -0.2650983643125788D4
      xinv(7,8) = -0.6610569847432101D3
      xinv(8,1) = 0.905943007623547D2
      xinv(8,2) = 0.4532784844440925D3
      xinv(8,3) = 0.3977386507582565D3
      xinv(8,4) = -0.9416114359647037D3
      xinv(8,5) = 0.961358754901498D1
      xinv(8,6) = 0.1195770634536251D3
      xinv(8,7) = 0.2499373357191282D3
      xinv(8,8) = 0.6359831792143064D2

      goto 1000

 500  continue
      xinv(1,1) = -0.409475795253919D5
      xinv(1,2) = -0.1625119519153035D6
      xinv(1,3) = -0.1221601301261945D6
      xinv(1,4) = -0.6384655723444662D5
      xinv(1,5) = 0.3894672188013366D6
      xinv(1,6) = -0.292978596250113D4
      xinv(1,7) = -0.3453015656557736D5
      xinv(1,8) = -0.8403029987245931D5
      xinv(1,9) = -0.8079722097070688D5
      xinv(1,10) = -0.1597560269451346D5
      xinv(2,1) = 0.2447025622926167D6
      xinv(2,2) = 0.1015958518771659D7
      xinv(2,3) = 0.7998149226969885D6
      xinv(2,4) = 0.4298519070175431D6
      xinv(2,5) = -0.2490327910778807D7
      xinv(2,6) = 0.172987633070631D5
      xinv(2,7) = 0.2113682204120962D6
      xinv(2,8) = 0.5271195244011894D6
      xinv(2,9) = 0.5144262988238237D6
      xinv(2,10) = 0.102314886409861D6
      xinv(3,1) = -0.644148066341683D6
      xinv(3,2) = -0.2792965584350808D7
      xinv(3,3) = -0.2305257785008599D7
      xinv(3,4) = -0.1276530475935926D7
      xinv(3,5) = 0.7018901911637016D7
      xinv(3,6) = -0.4506048845393519D5
      xinv(3,7) = -0.5694200304952046D6
      xinv(3,8) = -0.1455685753053223D7
      xinv(3,9) = -0.1443197007728059D7
      xinv(3,10) = -0.2888764852839543D6
      xinv(4,1) = 0.9808781429986706D6
      xinv(4,2) = 0.4432404806876069D7
      xinv(4,3) = 0.3838315576272213D7
      xinv(4,4) = 0.2194356038934521D7
      xinv(4,5) = -0.1144595456508147D8
      xinv(4,6) = 0.6798554366660939D5
      xinv(4,7) = 0.8863896177208164D6
      xinv(4,8) = 0.2323084813458226D7
      xinv(4,9) = 0.2341755202366627D7
      xinv(4,10) = 0.4719817001333924D6
      xinv(5,1) = -0.9526620900457856D6
      xinv(5,2) = -0.4476286758588098D7
      xinv(5,3) = -0.40681179778208D7
      xinv(5,4) = -0.2405924939945288D7
      xinv(5,5) = 0.1190299176639997D8
      xinv(5,6) = -0.6549428602217757D5
      xinv(5,7) = -0.8789754982083279D6
      xinv(5,8) = -0.2361412249941344D7
      xinv(5,9) = -0.242224212675238D7
      xinv(5,10) = -0.4918509629983634D6
      xinv(6,1) = 0.6122682828689225D6
      xinv(6,2) = 0.2984327254245365D7
      xinv(6,3) = 0.2846034224748765D7
      xinv(6,4) = 0.1744589820838437D7
      xinv(6,5) = -0.818721958270149D7
      xinv(6,6) = 0.4178987704590344D5
      xinv(6,7) = 0.576042663059888D6
      xinv(6,8) = 0.1585881550729814D7
      xinv(6,9) = 0.1656554800913041D7
      xinv(6,10) = 0.3390760123047862D6
      xinv(7,1) = -0.2604921260736874D6
      xinv(7,2) = -0.1313998566753303D7
      xinv(7,3) = -0.1314245615661388D7
      xinv(7,4) = -0.8365457658055928D6
      xinv(7,5) = 0.372528207429397D7
      xinv(7,6) = -0.1766579354125239D5
      xinv(7,7) = -0.2495923943741827D6
      xinv(7,8) = -0.703823799974276D6
      xinv(7,9) = -0.7491541434239853D6
      xinv(7,10) = -0.1546632097020928D6
      xinv(8,1) = 0.7077146835826727D5
      xinv(8,2) = 0.3685946195961626D6
      xinv(8,3) = 0.3863090375487359D6
      xinv(8,4) = 0.2557599881911758D6
      xinv(8,5) = -0.1081435113694342D7
      xinv(8,6) = 0.4772026851795869D4
      xinv(8,7) = 0.689736678650375D5
      xinv(8,8) = 0.1990994915820471D6
      xinv(8,9) = 0.2160679708027538D6
      xinv(8,10) = 0.4501784610481079D5
      xinv(9,1) = -0.1114491595488314D5
      xinv(9,2) = -0.5979909670968722D5
      xinv(9,3) = -0.6559541101761499D5
      xinv(9,4) = -0.4523640567885096D5
      xinv(9,5) = 0.1817758293610363D6
      xinv(9,6) = -0.7476267686651723D3
      xinv(9,7) = -0.1103516087047956D5
      xinv(9,8) = -0.3258464218798142D5
      xinv(9,9) = -0.3606964318138988D5
      xinv(9,10) = -0.7588768171611588D4
      xinv(10,1) = 0.7753214229539828D3
      xinv(10,2) = 0.4276758827943425D4
      xinv(10,3) = 0.4903158367893109D4
      xinv(10,4) = 0.3526389618427347D4
      xinv(10,5) = -0.1348162823721786D5
      xinv(10,6) = 0.5176987715964296D2
      xinv(10,7) = 0.7790714559339035D3
      xinv(10,8) = 0.2351364858006837D4
      xinv(10,9) = 0.2655869150276276D4
      xinv(10,10) = 0.5645838976852591D3

      goto 1000

 1000 continue

      do i=1,s
c compute y(alpha):

         do ii=1,s
            do j=0,2*s-1
               y(ii,j+1)=(alpha*aa(i)*c(ii)+1d0)**j
            end do
         end do

         call dgemm('n','n',s,2*s,2*s,1d0,y,s,xinv,2*s,0d0,tmp,s)

         call dcopy(s*s,tmp(1,1),1,vss(1,1,i),1)
         call dcopy(s*s,tmp(1,s+1),1,wss(1,1,i),1)


      end do
      return

      end



