holger@0: function betaDivPerElem = betaDivergencePerElement(x, y, beta)
holger@0: 
holger@0: switch beta
holger@0:     case 0
holger@0:         betaDivPerElem = x./y - log(x./y) - 1;
holger@0: 
holger@0:     case 1
holger@0:         betaDivPerElem = x .* log(x./y) + x - y;
holger@0: 
holger@0:     otherwise
holger@0:         betaDivPerElem = (x.^beta) / (beta*(beta-1)) ...
holger@0:                          + (y.^beta) / beta ...
holger@0:                          - (x .* (y.^(beta-1))) / (beta-1);
holger@0: end