There are several errors in your code. Consider this (I put comments where corrections were needed)
n=0;
eps=1;
a=0.1; %You need a way smaller parameter to converge!
x=[1;1];
A = [6 -3/2 ; -3/2 8]; %You have a bilinear positive definite form,
%you may use matrix form for convenience
while eps>1e-12 && n<100 %You had wrong termination conditions!!
gradf=2*A*x; %(gradf in terms of matrix)
f=x'*A*x; %you need to update f every iteration!!
eps=(norm(gradf)/(1+abs(f)))
disp(eps > 1e-12)
x=x-a*gradf;
%Now you can see the orbit towards minimum
plot(x(1),x(2),'o'); hold on
n=n+1;
end
n
x
eps
for instance with the current value a=.1 I get
n = 100
eps = 1.2308e-011
x =
1.0e-012 *
-0.2509
0.4688
That is I had to perform 100 iteration because my epsilon is still above the threshold. If I allow 200 iterations I get
n = 110
eps =
7.9705e-013
x =
1.0e-013 *
-0.1625
0.3036
I.e. 110 iterations are sufficient.

Case of a general f (i.e. not a quadratic form).
You can use, for instance, function handles, i.e. you define (before the while)
foo = @(x) 6*x(1)^2+8*x(2)^2-3*x(1)*x(2);
foo_x = @(x) 12*x(1)-3*x(2);
foo_y = @(x) 16*x(2)-3*x(1);
then, in the while you substitute
gradf = [foo_x(x);foo_y(x)];
f = foo(x);
P.S. for what concerns the while cycle, please keep in mind that you keep on iterating while you are not satisfied with your precision (eps>1e-12) AND your total number of iteration is below a given threshold (n<100).
Consider also that you are working in finite precision: a numerical algorithm can never reach the analytic solution (i.e. what you have with infinite precision and infinite iterations), therefore, you always have to set a threshold (eps, which should be above the machine precision \approx1e-16) and that is your 0.