Why does adding L1 penalty to R's optim slows things down so much (relative to no penalty or L2)? I'm running some optimizations with optim's implementation of BFGS.  The objective function is actually a computational algorithm, not just math.   I found when I add an L1 penalty, things slow down quite a bit.  Why might this be?  Is there something about L1 that slows things down?  Then how is it glmnet implementation of LASSO so fast?  
A quick Google search turned up a package call "lbfgs" which "finds the optimum of an objective plus the L1 norm of the problem’s parameters" with "a fast and memory-efficient implementation of these optimization routines, which is particularly suited for high-dimensional problems."  Should I be looking for solutions like this?
 A: I would guess that the reason why adding an L1 penalty slows things down significantly is that an L1 penalty is not differentiable (i.e. absolute value), while the L2 penalty is. This means that the surface of the function will not be smooth, and so standard quasi-Newton's methods will have a lot of trouble with this problems. Recall that one way to think of a quasi-Newton's method is that it makes a quadratic approximation of the function and then the initial proposal will the maximum of that approximation. If the quadratic approximation matches fairly well to the target function, we should expect the proposal to be close the maximum (or minimum, depending on how you look at the world). But if your target function is non-differentialable, this quadratic approximation may be very bad, thus taking many more iterations to converge (assuming it does). 
If you've found an R-package that implements BFGS for L1 penalties, by all means try it. BFGS, in general, is a very generic algorithm for optimization. As is the case with any generic algorithm, there will be plenty of special cases where it does not do well. Algorithms that are specially tailored to your problem clearly should do better (assuming the package is as good as it's author claims: I haven't heard of lbfgs, but there's a whole lot of great things I haven't heard of. Update: I have used R's lbfgs package, and the L-BFGS implementation they have is quite good! I still haven't used it's OWL-QN algorithm, which is what the OP is referring to). 
If it doesn't work out for you, you might want to try the "Nelder-Mead" method with R's optim. It does not use derivatives for optimization. As such, it will typically be slower on a smooth function but stabler on an unsmooth function such as you have.
A: I don't know why your problem slows down when you add an $L_1$ penalty. It probably depends on (1) what the problem is; (2) how you've coded it; and (3) the optimization method you're using.
I think the there's an "unspoken answer" to your question: the most efficient solutions to numerical problems are often tailor-made. General-purpose algorithms are just that: general purpose. Specialized solutions to specific problems will tend to work better, because we can bring to bear observations about how that particular problem is presented and its specific properties which are known to the analyst. To your specific question about glmnet, it has a number "tricks" which makes it highly efficient - for the particular problem that it's trying to solve! The Journal of Statistical Software paper on provides details: 


*

*Its optimization for all models (elastic net, ridge regression and not just LASSO) uses cyclical coordinate descent, which is a pretty good way to go about solving this problem. 

*The coefficients are computed along paths for a range of $\lambda$ values. So rather than wandering over the response surface for a single value of the regularization parameter $\lambda$, the it moves from largest to smallest values, using coefficient estimates from previous solutions as starting points. This exploits the fact that coefficient estimates ascend from smaller to larger values as $\lambda$ decreases; it doesn't have to re-solve the same problem over and over again from randomly-initialized starts as one would with a naive implementation of a standard optimization routine.


And it's coded in FORTRAN.
L-BFGS is a limited memory BFGS algorithm. While it has tricks that can make it more efficient than standard BFGS for some problems, it's not clear whether the problems that it solves have any bearing on the particular problem at hand. L-BFGS is one of the options in optim as well, so I'm not sure why you'd need an additional package.
Note that BFGS depends on derivatives, which are computed by finite differences when analytical forms are not provided. This could be where you get problems, because the $L_1$ penalty is not differentiable everywhere. Not only does this mean that you're probably not going to estimate LASSO coefficients at precisely 0, it might wreak havoc with updating from iteration to iteration.
