How to drive shape keys based on the angle between two bones?
For an ultra-fast driver for corrective shape keys based on the angle between two bones, you can use this one:



The driver type is Scripted Expression, with this line (depending on your Blender version):

Blender 2.7x:

poseBones[A].matrix.col[1] * poseBones[B].matrix.col[1]


Blender 2.8x:

(Same line as above, but replace * with @)


(In order for this expression to work you need to enable the Auto Run Python Scripts setting in Blender. On 2.80 it's in Preferences -> Save &amp; Load -> Auto Run Python Scripts and on 2.79 it's in User Preferences -> File -> Auto Run Python Scripts.)  

This expression uses the dot product of the bones' aim vectors, a mechanism known as a "pose reader", used with pose space deformation (PSD).
As the bones rotate during animation or posing, the driver will output a value in the range [-1, +1], going from -1 when the bones point in opposite directions, up to 0 when they become perpendicular, up to +1 when they point in the same direction.

Replace A and B in that expression with the indices of the two bones involved. The indices are from the bones list of the armature, and the only way to find them is by using Python. But it's simple:  


Select your armature object;
Change the editor to Python Console;
Type the next line in the console and hit Enter:


C.object.pose.bones.find("My Bone")


Instead of "My Bone", use the name of the bone whose index you want to find, then press Enter to see the result. Do this once for each bone so you'll have two numbers, like 12 and 103 for example. Put them in the expression (in any order), inside the poseBones brackets:

poseBones[103].matrix.col[1] * poseBones[12].matrix.col[1]


Next, the other driver variables are:


A Single Property type variable named 'poseBones', used in that expression. The ID type is set to Object, and the object itself is the armature that has the two bones involved. The actual path to the property in the object is pose.bones, write this in that field. See more details in the Notes section below.
(Blender 2.79 only, not needed for 2.80+):
A dummy Transform Channel type variable, it needs to point to any bone in that armature (any indeed, like the root bone). This variable is unused but is an essential piece, it's a way to force Blender 2.79 to evaluate the driver every time a bone in your armature moves from animation or posing.  Set the channel type and space to Location X and World Space, respectively, which makes this very fast to evaluate.
Again, this variable is not used in the expression and is only needed when on Blender 2.7x.


Your driver is ready to be used, with no roll changes or sudden flips.





Notes:


The two bones used in the expression should not be animated in scale (be it their own or inherited scale from their parents), they should only be animated in position and rotation. This means that they need to have a scale of exactly 1 in all axes, for the entire animation, otherwise this will mess with the speed of the driver because it affects the length of the matrix.col[1] vector.  


If you need to animate them in scale, simply add a child bone to them and animate the scale of the child and deform the mesh with that child bone, not the parent.  

You can adjust the timing and speed of this driver by playing with the driver graph curve, like using easing or bezier handles.


For example, when used with shape keys, you can make the graph minimum be zero and maximum be 1 (click to enlarge):

The natural rate of change of this driver is senoidal because of the dot product that it uses in the expression. If you happen to need linear speed, adjust the bézier handles of the curve in this way to keep the rate of change constant (click to enlarge):


The poseBones variable in the driver is referencing the .pose.bones property of the armature object (like bpy.data.objects("Armature")) and not the armature data (like bpy.data.armatures("Armature")). That .pose property exists in objects only.