implicit none
      real*8 t0,t1,t2,g,p(10000),s,pieq,a0,a1
      real*8 h,r,pk,x1,x2,rh,ss(0:10000)
      real*8 q(10000),acc,res(0:10000),sh
      integer i,k,n,l0,lp,le,l0h,lph,leh,k1,k2,l,lma
      
c      parameter(t0=0.0009,t1=0.000044,t2=0.090671,l0=30,lp=65,
c     1le=90,n=120,s=50,h=0.25,r=0.02)
c      parameter(t0=0.000204,t1=0.0000068,t2=0.111012,l0=35,lp=65,
c     1le=120,n=120,s=1,h=1,r=0.03)
      parameter(t0=0.000309,t1=0.0000231,t2=0.10005,l0=35,lp=65,
     1le=120,n=120,s=1,h=0.0833333,r=0.03)

      open(unit=10,file='acount.res')
      open(unit=15,file='reserve.res')

c-----------------------------------------

c        Aldre, rente og pensjon (?rlige) konverteres til h-enheter der h er 
c        h er tidsinkrementet     

      lph=lp/h
      leh=le/h
      l0h=l0/h
      rh=dexp(dlog(1.0+r)*h)-1
      sh=h*s

c       Gomperz-makeham overlevelsessannsynligheter beregnes
c       Funskjonen g er en hjelperutine, se bakerst p? denne progamfilen 

      do i=1,leh
        x1=(i-1)*h
        x2=i*h
        p(i)=dexp(-g(x2,t0,t1,t2)+g(x1,t0,t1,t2))
      enddo
      
c-------------------------------------------

c        Ekvivalenspremiebereregninger

      a0=0
      pk=1
      do k=0,lph-l0h-1
        a0=a0+pk
        pk=pk*p(l0h+k+1)/(1+rh)
      enddo
      a1=0
      do k=lph-l0h,leh-l0h-1
        a1=a1+pk
        pk=pk*p(l0h+k+1)/(1+rh)
      enddo
      pieq=sh*a1/a0
      write(*,100)a0,a1,pieq

c       Husk at pieq/h svarer til  en ?rlig premie. Kj?res programmer med h=1
c       er det den ?rlige ekvivalenspremien som bestemmes 

 100  format(1x,10f12.3)
 101  format(1x,i7,2f14.4)
c-------------------------------------

c      Laster inn i betalingsfunksjon

      do k=0,lph-l0h-1
         ss(k)=-pieq
      enddo
      do k=lph-l0h,leh-l0h
        ss(k)=sh
      enddo

c        her beregnes reserve p? tidspunkt k

      do k=0,leh-l0h
        res(k)=0
        pk=1
        do i=0,leh-l0h-k-1
          res(k)=res(k)+pk*ss(i+k)
          pk=pk*p(l0h+k+1+i)/(1.0+rh)
        enddo
        write(15,100)k*h,res(k)
      enddo
        
c----------------------------------------------

c       Her beregnes akkumulert spareverdie av det som er betalt inn 
c       p? tidspunkt k

      acc=-ss(0)
      k=0
      write(10,101)k,acc     
      do k=1,leh-l0h-1
           acc=acc*(1+rh)-ss(k)
        write(10,100)k*h,acc
      enddo

      end
c-------------------------------
      real*8 function g(x,t0,t1,t2)
      implicit none
      real*8 x,t0,t1,t2
     
      g=t0*x+(t1/t2)*(dexp(t2*x)-1)

      end