/*  solve differential equation for continuous time search with limit T on search time 
    agents still live forever. If no job accepted by time T then individual recieves
    b/rho */
    
 closeall;
 
 output file = search-de.out reset;
 
 /* set parameter values */
 
 rho = .05;
 b = 1;
 lambda = .3;
 t = 200;
 
 m = 10000;  @ number of intervals @
 
 gwstar = b;
 
 _intord = 40;
 
 proc ifp(y);
 retp( (y - gwstar) .* ( .2 * exp( - .2 * y) ) );
 endp;
 
 
 delta = t/(m-1);
 vvec = zeros(m,1);
 timeseq = seqa(0,delta,m);
 
 proc update(x);
 gwstar = rho*x;
 retp(  b - rho*x + lambda/rho * intquad1(&ifp,( 100 | (rho*x) )) ) ;
 endp;
 
 vvec[1] = b/rho;
 i = 2;
 do until i gt m;
   vvec[i] = vvec[i-1] + delta*update(vvec[i-1]);
   i = i+1;
 endo;
 
 library pgraph;
 graphset;
 @xy(timeseq,vvec*rho);@
 
 
 /** we can extend this to replicate the model van den berg looked at,
     namely, let individual receive a unemployment utility flow of 1 after his
     benefits have been exhausted
     assume that with benefits b = 3, and that unemployment benefits are lost after 9 periods */
     b = 3;
 /* now from the analysis above, we already have the constant steady state reservation wage after
    benefits have been exhausted */
    
    wstar_ss = rho*vvec[10000];
    gwstar = wstar_ss;
    t = 9;   /* maximum amount of time with unemployment benefits */
    m = 1000;
    delta = t/(m-1);
    vvec = zeros(m,1);
    vvec[1] = gwstar/rho;
     i = 2;
 do until i gt m;
   vvec[i] = vvec[i-1] + delta*update(vvec[i-1]);
   i = i+1;
 endo;   
 
 nvec = rev(rho*vvec) | (wstar_ss*ones(m,1));
 tseq = seqa(0,delta,2*m);
 
 graphset;
 xy(tseq,nvec);
 
 
 end;