fmincon is not applicable since you have integrality constraints.
The model can easily be converted to a mixed-integer linear problem (if c is non-negative using a standard epigraph reformulation of the max operator, otherwise using a big-M model of it), for which there are many solvers available for MATLAB
The MATLAB toolbox YALMIP makes the modelling easy. Install a MILP solver for good performance (otherwise it uses its own naive implementation)
N = 10;
c = rand(N,1);
p = randn(N,1);
a = intvar(N,1);
solvesdp([-5 <= a <= 5],p'*a+c'*max(a-1,0))
double(a)
Note though, I think you can solve this problem analytically.
