fmincon to solve constrained optimization problem with a nonlinear constraint

I've got a challenging problem to solve using the fmincon function. I've attached the question and posted the code from the question I solved prior to this one to give you a better idea of what I'm doing. Thanks for your guidance...here's the previous code.

% Used for plotting as it contains both x and y
z2 = @(x,y) (x/(sqrt(x.^2 + y.^2))).*(besselj(1, 3.8316.*sqrt(x.^2 + y.^2)));

% Used for finding min at each function
z0 = @(x) (x(1)/(sqrt(x(1).^2 + x(2).^2))).*(besselj(1, 3.8316.*sqrt(x(1).^2 + x(2).^2)));

% Used to find max at each function due to multiplying by (-1)
z1 = @(x) (-1)*(x(1)/(sqrt(x(1).^2 + x(2).^2))).*(besselj(1, 3.8316.*sqrt(x(1).^2 + x(2).^2)));

% Defining initial conditions
a = [.5, -.5];

% Running min/max search using fminsearch function
[xi,fvali,exitflagi,outputi] = fminsearch(z0, a);
[xa,fvala,exitflaga,outputa] = fminsearch(z1, a);

% Running min/max search using fminunc function
[xi1,fvali1,exitflagi1,outputi1] = fminunc(z0, a);
[xa1,fvala1,exitflaga1,outputa1] = fminunc(z1, a);

% Plotting a max found from the search
hold on
ezmeshc(z2, [-1, 1, -1, 1]);
plot3(xa(1,1),xa(1,2), fvala  , 'ob', 'MarkerSize', 12);
axis([-1 1 -1 1 -1 1]); % Expanded the axis to show marker
xlabel('x-axis');
ylabel('y-axis');
zlabel('z-axis');
legend('grid', 'rings', 'marker');
view(155, 9); % Adjusting starting orientation of plot

[BlkHoleSun]

pressed on with a solution, but not confident I'm correct. All the parts are on the plot, but the max I found using the fmincon function doesn't seem right. below is the separate function code and problem code I've written for the problem:

Function Code:

function [c,ceq] = nonlcon(x)
%Non-linearinequalityconstraints:
c = (0.6)^2-(x(1)-0.4)^2 + (x(2)-0.4)^2; 
ceq = [];

Problem Code:

% Used for plotting as it contains both x and y
z2 = @(x,y) (x/(sqrt(x.^2 + y.^2))).*(besselj(1, 3.8316.*sqrt(x.^2 + y.^2)));

% Used to find max at each function due to multiplying by (-1)
z1 = @(x) (-1)*(x(1)/(sqrt(x(1).^2 + x(2).^2))).*(besselj(1, 3.8316.*sqrt(x(1).^2 + x(2).^2)));
x0 = [0.5,-0.5];

[x,fval,exitflag,output] = fmincon(z1,x0,[],[],[],[],[],[],@nonlcon)

hold on
plot3(x(1,1),x(1,2), fval  , 'ob', 'MarkerSize', 12);
ezmeshc(z2, [-1, 1, -1, 1]);
[x,y,z] = cylinder(0.6, 40);
mesh(x+0.4, y+0.4, 1.2*z-0.6, 'FaceAlpha', 0.8)
axis equal
view(-64, 14);
xlabel('x-axis');
ylabel('y-axis');
zlabel('z-axis');
legend('marker', 'grid', 'rings');

Answer 1

[BlkHoleSun]

%%Set up inline functions to be used for fmincon and plotting
clear all
% Used for plotting as it contains both x and y
z2 = @(x,y) (x/(sqrt(x^2 + y^2)))*(besselj(1, 3.8316.*sqrt(x^2 + y^2)));

% Used to find max due to multiplying by (-1)
z0 = @(x) ((x(1))/(sqrt((x(1))^2 + (x(2))^2)))*(besselj(1, 3.8316*(sqrt((x(1))^2 + (x(2))^2))));
z1 = @(x) -1*(x(1)/(sqrt((x(1)^2) + (x(2)^2))))*(besselj(1, (3.8316*(sqrt((x(1)^2)+(x(2)^2))))));

%%Using fmincon function with embedded nonlcon function to find max
a = [.5 ,-.5]; %Starting point from problem 3
[xa,fvala,exitflaga,outputa] = fmincon(z1,a,[],[],[],[],[],[],@nonlcon,[]);
maxit = outputa.iterations;
fprintf('The max number of iterations using fminunc is %2.0f. \n', maxit)
%%Creating plot with problem 3 function, constraint and max found
hold on
plot3(xa(1),xa(2), -fvala , 'ob', 'MarkerSize', 12);
ezmeshc(z2, [-1, 1, -1, 1]);
[x,y,z] = cylinder(0.6, 40);
mesh(x+0.4, y+0.4, 1.2*z-0.6, 'FaceAlpha', 0.8)
axis equal
view(-64, 14);
xlabel('x-axis'); ylabel('y-axis'); zlabel('z-axis');
title('z(x,y) vs max found vs cylinder constraint')
legend('marker', 'grid', 'rings');