anchor.diff:          
   
       write(*,*) "Starting Solver"
  -    call system_clock(startTime, clock_rate, clock_max)
  +    call system_clock(COUNT_RATE=clock_rate, COUNT_MAX=clock_max)
  +    call system_clock(startTime)
       call solveOrder2(tEnd,nOuts,tDel,n,row,col,M,C,yLen,xDel,yDel,&
          c0,m0,y,nu,gam,alpha,beta,k,delta,ndiag,e1,e2,nit,esoln,dSelect,&
          fSelect,gSelect,avgIters,maxIters,avgNit,maxNit)
  -    call system_clock(endTime, clock_rate, clock_max)
  -    
  +    call system_clock(endTime)
  +   
  +    ! Convert times into seconds
  +    elapsedTime = real(endTime - startTime)/real(clock_rate)
  +
  +
   !    call solveOrder1(tEnd,nOuts,tDel,n,row,col,M,C,xLen,yLen,xDel,yDel&
   !       ,c0,m0,nu,gam,alpha,beta,k,delta,ndiag,e1,e2,nit,dSelect,fSelect)
     
  
positive.code:     !!$omp parallel do shared(MatrixM,yCof,diff,Mrhs,tDel,M,C,k,yConst,nu,m0,c0,gam,fSelect) private(i, f)
      !$omp parallel do private(f)
      do i = 1,n
          if (i .GE. n-row+1) then
              MatrixM(i,1) = MatrixM(i,1) - yCof*0.5*(diff(i-col)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + yCof*0.5*(diff(i-col)+diff(i))
          else
              MatrixM(i,5) = MatrixM(i,5) - yCof*0.5*(diff(i+col)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + yCof*0.5*(diff(i+col)+diff(i))
          endif
            
          if (MOD(i,col) == 0) then
              MatrixM(i,2) = MatrixM(i,2) - xCof*0.5*(diff(i-1)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + xCof*0.5*(diff(i-1)+diff(i))
          else
              MatrixM(i,4) = MatrixM(i,4) - xCof*0.5*(diff(i+1)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + xCof*0.5*(diff(i+1)+diff(i))
          endif
  
          if (MOD(i,col) == 1) then
              MatrixM(i,4) = MatrixM(i,4) - xCof*0.5*(diff(i+1)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + xCof*0.5*(diff(i+1)+diff(i))
          else
              MatrixM(i,2) = MatrixM(i,2) - xCof*0.5*(diff(i-1)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + xCof*0.5*(diff(i-1)+diff(i))
          endif
            
          if (i .LE. row) then
              MatrixM(i,5) = MatrixM(i,5) - yCof*0.5*(diff(i+col)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + yCof*0.5*(diff(i+col)+diff(i))
          else
              MatrixM(i,1) = MatrixM(i,1) - yCof*0.5*(diff(i-col)+diff(i))
              MatrixM(i,3) = MatrixM(i,3) + yCof*0.5*(diff(i-col)+diff(i))
          endif
          
          call fFunc(M(i),C(i),f,k,yConst,nu,m0,c0,gam,fSelect)
  
negative.code:     MatrixM(:,:) = 0
  
      ! Compute all the diffusion coefficients
      !$omp parallel do shared(M,diff,delta,alpha,beta,dSelect)
      do i = 1,n
          call dFunc(M(i), diff(i), delta, alpha, beta, dSelect)
      enddo
  