***********************************************************************
      SUBROUTINE GAUSS_32(X,TS,TF,STEP,ERR,N,NOR,NI,NS,NBS,NF,FUN)
      ! Gauss-Everhart Method (Radau & Lobatto Spacings)
      ! Written by Avdyushev V. (2010)
      IMPLICIT REAL*8(A-H,O-Z)
      PARAMETER (IPS0=100,EPS=1.D-15,SGM=3.16D0)
      REAL*8, ALLOCATABLE  :: A(:,:),B(:,:),C(:,:),D(:,:),E(:,:)
      REAL*8, ALLOCATABLE  :: SMS(:,:),SPC(:),BL(:,:)
      REAL*8, DIMENSION(N) :: X,Y,F,F0,P,XK
      LOGICAL LS,VS,INC,CNV,RADAU
      COMMON/SPACING/SPACING(256) 
      EXTERNAL FUN

      IF(TS.EQ.TF) RETURN
      IF(NOR.LT.2.OR.NOR.GT.32) 
     *STOP 'Inadmissible Order!'

      K=NOR/2 
      ALLOCATE(SPC(0:K),SMS(K,K))
      ALLOCATE(A(K,N),B(K,N),BL(K,N))
      ALLOCATE(C(0:K,0:K),D(0:K,0:K),E(0:K,0:K))
      SPC(0)=0.D0; L=K*(NOR-K-1); SPC(1:K)=SPACING(L+1:L+K)

      DO I=1,K; DO J=1,I-1
      SMS(I,J)=1.D0/(SPC(I)-SPC(J))
      END DO; END DO

      C=0.D0; C(0,0)=1.D0
      D=0.D0; D(0,0)=1.D0
      DO J=1,K; DO I=J,K
      C(I,J)=C(I-1,J-1)-SPC(I-1)*C(I-1,J)
      D(I,J)=D(I-1,J-1)+SPC(J  )*D(I-1,J)
      END DO; END DO

      E=0.D0; E(0:K,0)=1.D0
      DO I=1,K; DO J=1,I
      E(I,J)=E(I-1,J-1)+E(I-1,J)
      END DO; END DO

      DO I=1,K
      FLT=DFLOAT(I+1)
      C(:,I)=C(:,I)/FLT
      D(I,:)=D(I,:)*FLT
      E(I,:)=E(I,:)*FLT
      END DO

      RADAU=NOR-K*2.EQ.1
      DEG=1.D0/DFLOAT(K+1)
      RMX=SGM**DEG; RMN=1.D0/RMX
      DIR=(TF-TS)/DABS(TF-TS); H_=DABS(STEP)
      VS=ERR.NE.0.D0; CNV=NI.LE.0; INC=.TRUE. 
      T=TS; NS=0; NBS=0; NST=0; NF=0; IPS=IPS0
    3 H=H_*DIR; R=1.D0; LS=.FALSE.
      B=0.D0; XK=0.D0

      IF(H_.EQ.0.D0) THEN
      IF(.NOT.VS)
     *STOP 'Constant Step is Zero!'
      CALL FUN(TS,X,F0); NF=NF+1
      H_=EPS; B_=0.D0
      DO WHILE(B_.EQ.0.D0)
      H_=H_*10.D0; H=H_*DIR; Y=X+H*F0
      IF(DABS(TF-TS).LE.H_) GOTO 1
      CALL FUN(TS+H,Y,F); NF=NF+1
      DO J=1,N; B_=B_+(F(J)-F0(J))**2; END DO 
      END DO
      H_=DSQRT(2.D0*ERR*H_/DSQRT(B_))
      H=H_*DIR
      END IF

    1 IF(DABS(TF-T).LE.H_) THEN
      STEP=H_; H=TF-T; H_=DABS(H) 
      R=R*H_/STEP; LS=.TRUE.; END IF

      IF(NS.LE.1) BL=B

      Q=1.D0
      DO L=1,K
      P=0.D0
      DO M=L,K
      P=P+E(M,L)*B(M,:)
      END DO 
      B(L,:)=B(L,:)-BL(L,:)
      Q=Q*R; BL(L,:)=P*Q/DFLOAT(L+1)
      B(L,:)=BL(L,:)+B(L,:)
      END DO

      DO L=1,K; P=0.D0; DO M=L,K
      P=P+D(M,L)*B(M,:)
      END DO; A(L,:)=P; END DO

      IF(RADAU.OR.NS.EQ.0) THEN 
      CALL FUN(T,X,F0); NF=NF+1
      ELSE; F0=F; END IF

      DO IT=1,IPS

      DO I=1,K

      S=SPC(I); P=0.D0
      DO J=K,1,-1; P=S*(B(J,:)+P); END DO
      Y=X+H*S*(F0+P)

      CALL FUN(T+H*S,Y,F); NF=NF+1

      P=(F-F0)/S
      DO J=1,I-1; P=(P-A(J,:))*SMS(I,J); END DO
      DO J=1,I; B(J,:)=B(J,:)+C(I,J)*(P-A(I,:)); END DO
      A(I,:)=P

      END DO

      IF(ALL(DABS(Y-XK).LE.DABS(EPS*Y))) EXIT
      XK=Y

      END DO

      IF(IT.GT.IPS0) THEN
      IF(CNV.AND.NS.NE.0) NBS=NBS+1
      END IF

      IF(VS) THEN
      B_=0.D0 
      DO J=1,N
      B_=B_+B(K,J)**2
      END DO 
      IF(B_.NE.0.D0) THEN
      R=(ERR/DSQRT(B_)/H_)**DEG
      ELSE; R=RMX; END IF
      IF(NS.EQ.0) THEN 
      INC=.NOT.LS.AND.INC
      IF(INC.AND.R.GT.RMX.OR.R.LT.RMN) THEN
      NST=NST+1; IF(NST.GT.10.AND.R.GT.RMX) GOTO 4
      H_=H_*R; GOTO 3 
      END IF; END IF
    4 IF(R.GT.RMX) R=RMX
      END IF

      IF(RADAU) THEN; P=0.D0
      DO J=K,1,-1; P=B(J,:)+P; END DO
      X=X+H*(F0+P); ELSE; X=Y; END IF
      T=T+H; NS=NS+1; H=H*R; H_=DABS(H) 
      XK=0.D0

      IF(CNV) THEN; IPS=IPS0
      ELSE; IPS=NI; END IF

      IF(LS) GOTO 2; GOTO 1

    2 DEALLOCATE(A,B,C,D,E,BL,SMS,SPC)

      END SUBROUTINE GAUSS_32
***********************************************************************
      BLOCK DATA
      REAL*8, DIMENSION(256) :: SPACING
      COMMON/SPACING/SPACING
      DATA SPACING/
     *1.00000000000000000000000000000000,
     *0.66666666666666666666666666666667,
     *0.50000000000000000000000000000000,
     *1.00000000000000000000000000000000,
     *0.35505102572168219018027159252941,
     *0.84494897427831780981972840747059,
     *0.27639320225002103035908263312687,
     *0.72360679774997896964091736687313,
     *1.00000000000000000000000000000000,
     *0.21234053823915294397475811012400,
     *0.59053313555926528913507374793117,
     *0.91141204048729605260445385623054,
     *0.17267316464601142810085377187657,
     *0.50000000000000000000000000000000,
     *0.82732683535398857189914622812343,
     *1.00000000000000000000000000000000,
     *0.13975986434378055215208708112488,
     *0.41640956763108317994330233133708,
     *0.72315698636187617231995400231437,
     *0.94289580388548231780687880744588,
     *0.11747233803526765357449851302033,
     *0.35738424175967745184292450297956,
     *0.64261575824032254815707549702044,
     *0.88252766196473234642550148697967,
     *1.00000000000000000000000000000000,
     *0.09853508579882642612349889788775,
     *0.30453572664636390548538517627883,
     *0.56202518975261385599498747999477,
     *0.80198658212639182746420786320470,
     *0.96019014294853125765919330990667,
     *0.08488805186071653506398389301626,
     *0.26557560326464289309811405904562,
     *0.50000000000000000000000000000000,
     *0.73442439673535710690188594095438,
     *0.91511194813928346493601610698373,
     *1.00000000000000000000000000000000,
     *0.07305432868025888514812603418031,
     *0.23076613796994549908311663988435,
     *0.44132848122844986791860665819448,
     *0.66301530971884570090294702791922,
     *0.85192140033151570815002314750402,
     *0.97068357284021510802794972308684,
     *0.06412992574519669233127711938966,
     *0.20414990928342884892774463430102,
     *0.39535039104876056561567136982732,
     *0.60464960895123943438432863017268,
     *0.79585009071657115107225536569898,
     *0.93587007425480330766872288061033,
     *1.00000000000000000000000000000000,
     *0.05626256053692214646565219103231,
     *0.18024069173689236498757994280918,
     *0.35262471711316963737390777017124,
     *0.54715362633055538300144855765235,
     *0.73421017721541053152321060830661,
     *0.88532094683909576809035976293249,
     *0.97752061356128750189117450042915,
     *0.05012100229426992134382737779083,
     *0.16140686024463112327705728645432,
     *0.31844126808691092064462396564567,
     *0.50000000000000000000000000000000,
     *0.68155873191308907935537603435433,
     *0.83859313975536887672294271354567,
     *0.94987899770573007865617262220917,
     *1.00000000000000000000000000000000,
     *0.04463395528996985073312102185830,
     *0.14436625704214557148521852022821,
     *0.28682475714443051894868623974909,
     *0.45481331519657335096772777004679,
     *0.62806783541672769756914603951737,
     *0.78569152060436924164245873241833,
     *0.90867639210020604399625854192546,
     *0.98222008485263654818679489896232,
     *0.04023304591677059308553366958883,
     *0.13061306744724746249844691257008,
     *0.26103752509477775216941245363437,
     *0.41736052116680648768689011702091,
     *0.58263947883319351231310988297909,
     *0.73896247490522224783058754636563,
     *0.86938693255275253750155308742992,
     *0.95976695408322940691446633041117,
     *1.00000000000000000000000000000000,
     *0.03625781288320946094116430076808,
     *0.11807897878999870019228511199474,
     *0.23717698481496038531730669285327,
     *0.38188276530470597536077024839650,
     *0.53802959891898906511685689131944,
     *0.69033242007236218294037953277052,
     *0.82388334383700471813682425392743,
     *0.92561261029080395536408181404400,
     *0.98558759035112345136717325918919,
     *0.03299928479597043283386293195030,
     *0.10775826316842779068879109194577,
     *0.21738233650189749676451801526112,
     *0.35212093220653030428404424222047,
     *0.50000000000000000000000000000000,
     *0.64787906779346969571595575777953,
     *0.78261766349810250323548198473888,
     *0.89224173683157220931120890805423,
     *0.96700071520402956716613706804969,
     *1.00000000000000000000000000000000,
     *0.03002903216148649704306435763440,
     *0.09828901220985322965120102159023,
     *0.19902107896310115486205369838276,
     *0.32405553832333489264284949106524,
     *0.46326123428433936712690482228811,
     *0.60536015311421315703804789492239,
     *0.73884032399154375973394834194515,
     *0.85288855035692975957240056442018,
     *0.93826792812285187447737063280575,
     *0.98808238656758440309025441304101,
     *0.02755036388855888829620993084839,
     *0.09036033917799666082567920914154,
     *0.18356192348406966116879757277817,
     *0.30023452951732553386782510421652,
     *0.43172353357253622256796907213015,
     *0.56827646642746377743203092786985,
     *0.69976547048267446613217489578348,
     *0.81643807651593033883120242722183,
     *0.90963966082200333917432079085845,
     *0.97244963611144111170379006915161,
     *1.00000000000000000000000000000000,
     *0.02527362039752034975333118646162,
     *0.08304161344740514670686537298197,
     *0.16917510037718142596943345609434,
     *0.27779671510903207443667869219539,
     *0.40150272023286081677227928632696,
     *0.53186238691041595791688961924225,
     *0.65999184208533481176639476610298,
     *0.77715939295616214449216854654264,
     *0.87538077485555692626470041273609,
     *0.94796454887281944741645730422704,
     *0.98998171953831959415697527013220,
     *0.02334507667891804405154726762227,
     *0.07682621767406384156703719645062,
     *0.15690576545912128696362048021682,
     *0.25854508945433189912653138318153,
     *0.37535653494688000371566314981289,
     *0.50000000000000000000000000000000,
     *0.62464346505311999628433685018711,
     *0.74145491054566810087346861681847,
     *0.84309423454087871303637951978318,
     *0.92317378232593615843296280354938,
     *0.97665492332108195594845273237772,
     *1.00000000000000000000000000000000,
     *0.02156206316585036090809308308300,
     *0.07105789873558898215118984486568,
     *0.14544745623506411920978808359438,
     *0.24040111047477294625739742595658,
     *0.35039934972274500723370827665702,
     *0.46904915068718232937106975710318,
     *0.58945491879854231751096564944067,
     *0.70461911573741977837708316562399,
     *0.80784894547014595900894274377139,
     *0.89314550911652334236589322956750,
     *0.95555353684459227697453320121445,
     *0.99146094501157258063133553912216,
     *0.02003247736636954932244991899228,
     *0.06609947308482637449988989854586,
     *0.13556570045433692970766379973955,
     *0.22468029853567647234168864707046,
     *0.32863799332864357747804829817916,
     *0.44183406555814806617061164513192,
     *0.55816593444185193382938835486808,
     *0.67136200667135642252195170182084,
     *0.77531970146432352765831135292954,
     *0.86443429954566307029233620026044,
     *0.93390052691517362550011010145413,
     *0.97996752263363045067755008100771,
     *1.00000000000000000000000000000000,
     *0.01861036501098785143971937784028,
     *0.06147554089926898760236661323470,
     *0.12630517869331058063228543286826,
     *0.20984297172656251444713666750062,
     *0.30789899828039834310295804831233,
     *0.41555603597865954449577915218909,
     *0.52741561399588227482490535732140,
     *0.63786860271776119959131870177268,
     *0.74137645929423748341020926717731,
     *0.83274898860844226850447752124064,
     *0.90740477530099736471710862456138,
     *0.96160186126032164962316747513586,
     *0.99263534897391067834930850158618,
     *0.01737703674808071360207430396519,
     *0.05745897788851185058729918425888,
     *0.11824015502409239964794076201185,
     *0.19687339726507714443823503068163,
     *0.28968097264316375953905153063071,
     *0.39232302231810288088716027686354,
     *0.50000000000000000000000000000000,
     *0.60767697768189711911283972313646,
     *0.71031902735683624046094846936929,
     *0.80312660273492285556176496931837,
     *0.88175984497590760035205923798815,
     *0.94254102211148814941270081574111,
     *0.98262296325191928639792569603480,
     *1.00000000000000000000000000000000,
     *0.01622476590139976171877199085899,
     *0.05369729993972461646659405657527,
     *0.11065719118048446030912833905351,
     *0.18461026055652535802692578138794,
     *0.27232354711073531456397327200967,
     *0.36996331162959604211551940586857,
     *0.47326213866012696794627305518712,
     *0.57770534269242974215959101248190,
     *0.67872825601106382559765960258696,
     *0.77191572935074200846535590138006,
     *0.85319513231878627007633983473920,
     *0.91901450031804481560754869226002,
     *0.96649859546798685996403607142943,
     *0.99358323920718154318917953590549,
     *0.01521597686489103352387863081627,
     *0.05039973345326395350268586924007,
     *0.10399585406909246803445586451842,
     *0.17380564855875345526605839017970,
     *0.25697028905643119410905460707656,
     *0.35008476554961839595082327263885,
     *0.44933686323902527607848349747704,
     *0.55066313676097472392151650252296,
     *0.64991523445038160404917672736115,
     *0.74302971094356880589094539292344,
     *0.82619435144124654473394160982029,
     *0.89600414593090753196554413548157,
     *0.94960026654673604649731413075992,
     *0.98478402313510896647612136918373,
     *1.00000000000000000000000000000000,
     *0.01426945473682577473409936694087,
     *0.04729959009416668566195579247573,
     *0.09771329932062197336876149533799,
     *0.16356903939438987602444091434582,
     *0.24233526096865728800292572225971,
     *0.33098480497004012346130436094686,
     *0.42611083909331411932854614476247,
     *0.52405769153676513942741100798415,
     *0.62106131135302196189347099085723,
     *0.71339391374247294001597395451560,
     *0.79750724494989595243178001167977,
     *0.87016897444640894402874546190571,
     *0.92858704688484115994521609825327,
     *0.97051770135205751336835901528200,
     *0.99435931102748829024249353342056,
     *0.01343391168429084292151024906313,
     *0.04456000204221320218809874680113,
     *0.09215187438911484644662472338123,
     *0.15448550968615764730254032131377,
     *0.22930730033494923043813329624797,
     *0.31391278321726147904638265963237,
     *0.40524401324084130584786849262344,
     *0.50000000000000000000000000000000,
     *0.59475598675915869415213150737656,
     *0.68608721678273852095361734036763,
     *0.77069269966505076956186670375203,
     *0.84551449031384235269745967868623,
     *0.90784812561088515355337527661876,
     *0.95543999795778679781190125319886,
     *0.98656608831570915707848975093686,
     *1.00000000000000000000000000000000/
      END BLOCK DATA
***********************************************************************
      SUBROUTINE FUN(T,X,F)
      IMPLICIT REAL*8(A-H,O-Z)
      DIMENSION X(4),F(4)
      R=DSQRT(X(1)**2+X(2)**2)
      F(1:2)= X(3:4)
      F(3:4)=-X(1:2)/R**3
      END SUBROUTINE FUN
