subroutine gauss(n,a,b,x,w)
  ! subroutine gauss for gaussian integration   76.08.12.
  ! the integral of f(x) from a to b is given by the sum
  ! over i from 1 to n of w(i)*f(x(i)).
  ! maximum n is 48
  implicit real*8(a-h,o-z)
  implicit integer(i-n)
  dimension x(n),w(n),e06(3),w06(3),e10(5),w10(5),e12(6),w12(6), &
       e16(8),w16(8),e20(10),w20(10),e24(12),w24(12), &
       e32(16),w32(16),e40(20),w40(20),e48(24),w48(24)
  data e06/.238619186083197d0,.661209386466265d0,.932469514203152d0/
  data w06/.467913934572691d0,.360761573048139d0,.171324492379170d0/
  data e10/.148874338981631d0,.433395394129247d0,.679409568299024d0, &
       .865063366688985d0,.973906528517172d0/
  data w10/.295524224714753d0,.269266719309996d0,.219086362515982d0, &
       .149451349150581d0,.066671344308688d0/
  data e12/.125233408511469d0,.367831498998180d0,.587317954286617d0, &
       .769902674194305d0,.904117256370475d0,.981560634246719d0/
  data w12/.249147045813403d0,.233492536538355d0,.203167426723066d0, &
       .160078328543346d0,.106939325995318d0,.047175336386512d0/
  data e16/.0950125098376374d0,.2816035507792589d0, &
       .4580167776572274d0,.6178762444026437d0, &
       .7554044083550030d0,.8656312023878317d0, &
       .9445750230732326d0,.9894009349916499d0/
  data w16/.1894506104550685d0,.1826034150449236d0, &
       .1691565193950025d0,.1495959888165767d0, &
       .1246289712555339d0,.0951585116824928d0, &
       .0622535239386479d0,.0271524594117541d0/
  data e20/.0765265211334973d0,.2277858511416451d0, &
       .3737060887154196d0,.5108670019508271d0, &
       .6360536807265150d0,.7463319064601508d0, &
       .8391169718222188d0,.9122344282513259d0, &
       .9639719272779138d0,.9931285991850949d0/
  data w20/.1527533871307259d0,.1491729864726037d0, &
       .1420961093183821d0,.1316886384491766d0, &
       .1181945319615184d0,.1019301198172404d0, &
       .0832767415767047d0,.0626720483341091d0, &
       .0406014298003869d0,.0176140071391521d0/
  data e24/.0640568928626056d0,.1911188674736163d0, &
       .3150426796961634d0,.4337935076260451d0, &
       .5454214713888395d0,.6480936519369756d0, &
       .7401241915785544d0,.8200019859739029d0, &
       .8864155270044010d0,.9382745520027328d0, &
       .9747285559713095d0,.9951872199970214d0/
  data w24/.1279381953467522d0,.1258374563468283d0, &
       .1216704729278034d0,.1155056680537256d0, &
       .1074442701159656d0,.0976186521041139d0, &
       .0861901615319533d0,.0733464814110803d0, &
       .0592985849154368d0,.0442774388174198d0, &
       .0285313886289337d0,.0123412297999872d0/
  data e32/.0483076656877383d0,.1444719615827965d0, &
       .2392873622521371d0,.3318686022821276d0, &
       .4213512761306353d0,.5068999089322294d0, &
       .5877157572407623d0,.6630442669302152d0, &
       .7321821187402897d0,.7944837959679424d0, &
       .8493676137325700d0,.8963211557660521d0, &
       .9349060759377397d0,.9647622555875064d0, &
       .9856115115452683d0,.9972638618494816d0/
  data w32/.0965400885147278d0,.0956387200792749d0, &
       .0938443990808046d0,.0911738786957639d0, &
       .0876520930044038d0,.0833119242269468d0, &
       .0781938957870703d0,.0723457941088485d0, &
       .0658222227763618d0,.0586840934785355d0, &
       .0509980592623762d0,.0428358980222267d0, &
       .0342738629130214d0,.0253920653092621d0, &
       .0162743947309057d0,.0070186100094701d0/
  data e40/.0387724175060508d0,.1160840706752552d0, &
       .1926975807013711d0,.2681521850072537d0, &
       .3419940908257585d0,.4137792043716050d0, &
       .4830758016861787d0,.5494671250951282d0, &
       .6125538896679802d0,.6719566846141795d0, &
       .7273182551899271d0,.7783056514265194d0, &
       .8246122308333117d0,.8659595032122595d0, &
       .9020988069688743d0,.9328128082786765d0, &
       .9579168192137917d0,.9772599499837743d0, &
       .9907262386994570d0,.9982377097105592d0/
  data w40/.0775059479784248d0,.0770398181642480d0, &
       .0761103619006262d0,.0747231690579683d0, &
       .0728865823958041d0,.0706116473912868d0, &
       .0679120458152339d0,.0648040134566010d0, &
       .0613062424929289d0,.0574397690993916d0, &
       .0532278469839368d0,.0486958076350722d0, &
       .0438709081856733d0,.0387821679744720d0, &
       .0334601952825478d0,.0279370069800234d0, &
       .0222458491941670d0,.0164210583819079d0, &
       .0104982845311528d0,.0045212770985332d0/
  data e48/ &
       0.032380170962869326033d0 , 0.097004699209462698930d0, &
       0.161222356068891718056d0 , 0.224763790394689061225d0, &
       0.287362487355455576736d0 , 0.348755886292160738160d0, &
       0.408686481990716729916d0 , 0.466902904750958404545d0, &
       0.523160974722233033678d0 , 0.577224726083972703818d0, &
       0.628867396776513623995d0 , 0.677872379632663905212d0, &
       0.724034130923814654674d0 , 0.767159032515740339254d0, &
       0.807066204029442627083d0 , 0.843588261624393530711d0, &
       0.876572020274247885906d0 , 0.905879136715569672822d0, &
       0.931386690706554333114d0 , 0.952987703160430860723d0, &
       0.970591592546247250461d0 , 0.984124583722826857745d0, &
       0.993530172266350757548d0 , 0.998771007252426118601d0/
  data w48/ &
       0.064737696812683922503d0 , 0.064466164435950082207d0, &
       0.063924238584648186624d0 , 0.063114192286254025567d0, &
       0.062039423159892663904d0 , 0.060704439165893880053d0, &
       0.059114839698398635746d0 , 0.057277292100403215705d0, &
       0.055199503699984162868d0 , 0.052890189485193667096d0, &
       0.050359035553854474958d0 , 0.047616658492490474826d0, &
       0.044674560856694280419d0 , 0.041545082943464749214d0, &
       0.038241351065830706317d0 , 0.034777222564770438893d0, &
       0.031167227832798088902d0 , 0.027426509708356948200d0, &
       0.023570760839324379141d0 , 0.019616160457355527814d0, &
       0.015579315722943848728d0 , 0.011477234579234539490d0, &
       0.007327553901276262102d0 , 0.003153346052305838633d0/
  e1 = .5d0*(b - a)
  e2 = e1 + a
  if(n == 6) go to  6
  if(n == 10) go to 10
  if(n == 12) go to 12
  if(n == 16) go to 16
  if(n == 20) go to 20
  if(n == 24) go to 24
  if(n == 32) go to 32
  if(n == 40) go to 40
  if(n == 48) go to 48
  xl = dsqrt(-e06(1))
  ! top2rx 2013.08.09 kino      stop 'gauss: wrong no. points'
  call rx( 'gauss: wrong no. points')
6 continue
  do  7 i = 1,3
     ak = e1*e06(i)
     x(i+3) = e2 + ak
     x(4-i) = e2 - ak
     w(i+3) = w06(i)*e1
     w(4-i) = w06(i)*e1
7 enddo
  return
10 continue
  do 11 i = 1,5
     ak = e1*e10(i)
     x(i+5) = e2 + ak
     x(6-i) = e2 - ak
     w(i+5) = w10(i)*e1
     w(6-i) = w10(i)*e1
11 enddo
  return
12 continue
  do 13 i = 1,6
     ak = e1*e12(i)
     x(i+6) = e2 + ak
     x(7-i) = e2 - ak
     w(i+6) = w12(i)*e1
     w(7-i) = w12(i)*e1
13 enddo
  return
16 continue
  do 17 i = 1,8
     ak = e1*e16(i)
     x(i+8) = e2 + ak
     x(9-i) = e2 - ak
     w(i+8) = w16(i)*e1
     w(9-i) = w16(i)*e1
17 enddo
  return
20 continue
  do 21 i = 1,10
     ak = e1*e20(i)
     x(i+10) = e2 + ak
     x(11-i) = e2 - ak
     w(i+10) = w20(i)*e1
     w(11-i) = w20(i)*e1
21 enddo
  return
24 continue
  do 25 i = 1,12
     ak = e1*e24(i)
     x(i+12) = e2 + ak
     x(13-i) = e2 - ak
     w(i+12) = w24(i)*e1
     w(13-i) = w24(i)*e1
25 enddo
  return
32 continue
  do 33 i = 1,16
     ak = e1*e32(i)
     x(i+16) = e2 + ak
     x(17-i) = e2 - ak
     w(i+16) = w32(i)*e1
     w(17-i) = w32(i)*e1
33 enddo
  return
40 continue
  do 41 i = 1,20
     ak = e1*e40(i)
     x(i+20) = e2 + ak
     x(21-i) = e2 - ak
     w(i+20) = w40(i)*e1
     w(21-i) = w40(i)*e1
41 enddo
  return
48 continue
  do 49 i = 1,24
     ak      = e1*e48(i)
     x(i+24) = e2 + ak
     x(25-i) = e2 - ak
     w(i+24) = w48(i)*e1
     w(25-i) = w48(i)*e1
49 enddo
  return
end subroutine gauss
