function [allslips,alltimes]=get_slip_hotfriction(prerate,costress,nafterpatches,stressG,stopt,sigAmB)

%Note: numsteps is number of time steps(to be determined by this algorithm)

% OUTPUTS:
%allslips is 2*nafterpatches x numsteps matrix of incremental slip during each
%   time step
%alltimes is "numsteps" length vector with time of each time step 

% INPUTS:
%pmpost is npatch x 7 matrix with model patch geometry where npatch is number
%   of patches
%prerate is 2*nafterpatches x 1 vector of pre-earthquake slip rates (1st npatch
%   values for strike-slip, second npatch values for dip-slip component)
%costress is 2*nafterpatches x 1 vector with coseismic shear stress change (1st npatch
%   values for strike-slip, second npatch values for dip-slip component
%nafterpatches is the number of slip patches
%stressG is 2*nafterpatches x 2*nafterpatches matrix the relates slip (s) to
%   shear stress (tau) on each patch, tau=stressG*s
%stopt is stop-time -- calculates slip from time=0 to time=stopt
%sigAmB is 2*nafterpatches x 1 vector of values sigma*(A-B) where A,B are
%   rate-state parameters and sigma is normal stress
 

%initialize some parameters
newstress=costress;
totalslip=0;
loop=0;
time=0;
delt=0;
alltimes=[];


%begin numerical integration by looping over time steps
while (time-delt)<stopt

    loop=loop+1;

    
    %sdot is slip rate
    sdot=prerate.*exp(newstress./[sigAmB;sigAmB]);
    
%abort if slip rate becoomes very high
  if max(sdot)>10^50
       totalslip=0*sdot; 
       sliprates=0*sdot;
       allslips=0*sdot;
       alltimes=0;
       disp('run-away')
      break
  end
  
   %delt is time step 
   %delt is function of sliprate
   delt=.01/max(sdot);
   
   %integrate overt time-step to get slip
   slip=sdot*delt;
   
   %stress due to slip during this time increment
   stress=stressG*slip;  
   
   
   %add incremement stress to total stress
   newstress=stress+newstress;
   
    
   %cumulative slip
   totalslip=totalslip+slip;
   
   %store sliprates, cumulative slip, times
   sliprates(:,loop)=sdot;
   allslips(:,loop)=slip;
   DELT(loop)=delt;
   alltimes(loop)=time;
   time=time+delt;

 
end