712.sagews

b, c = var('b,c')
L = plot3d(4*b^3*(2*b^2-3*b*c+b*c^2-c^2)-2*b*(4*b^3-c^3-3*b*c^2)*(1-c)*(b-c),
    (b,0,1), (c,0,1), adaptive=True, color=rainbow(60, 'rgbtuple'), max_depth=0)
Q = plot3d(0, (b,0,1), (c,0,1), color="lightblue", opacity=0.8)
L + Q