How can I fix this?

Whenever I run the program I get this:

??? In an assignment A(I) = B, the number of elements in B and I must be the same.

Error in ==> HW2P2ALMUSALLAM>funcalc at 100 f(1) = sin(x) + y.^2 + log(z) - 7;

Error in ==> HW2P2ALMUSALLAM>jacobcalc at 73 f = funcalc(x,y,z);

Error in ==> HW2P2ALMUSALLAM>newtonraphson at 36 J = jacobcalc(xo,yo,zo);

Error in ==> HW2P2ALMUSALLAM at 15 [x numit error] = newtonraphson(xo, yo, zo, rtol, maxit, rho);

Here is my program:

function HW2P2ALMUSALLAM

    % Solve exp(x) = sin(Pi*x/3) using Newton-Raphson method.

    % Specify parameters for Newton-Raphson method.
    rtol = 1e-5;
    maxit = 100;
    rho = 1;
    xo = 1;
    yo = 1;
    zo = 1;

    % Solve by Newton-Raphson.
    [x numit error] = newtonraphson(xo, yo, zo, rtol, maxit, rho);

    % Output results using fprintf (formatted output)
    fprintf('\n')
    fprintf('x = %8.5f \n',x)
    fprintf('Number of iterations = %3.0f \n',numit)
    fprintf('Relative error = %8.5e \n',error)

end

function [x numit error] = newtonraphson(xo, yo, zo, rtol, maxit, rho)

    % Initialize Iteration Variables
    error = 2*rtol;
    numit = 0;

    % Newton-Raphson Iteration Loop
    while error > rtol && numit < maxit

        % Compute Jacobian Matrix
        J = jacobcalc(xo,yo,zo);

        % Evaluate Function
        f = funcalc(xo,yo,zo);

        % Solve for delta x
        deltax = J\(-f); % J\ -> Inverse Jacobian or use inv(J)

        % Update x vector
        x = xo + rho * deltax;
        y = xo + rho * deltax;
        z = xo + rho * deltax;

        % Compute Error
        error = max( abs( (x - xo)./(x + eps) ) );

        % Increment loop counter & update xo value
        numit = numit + 1;
        xo = x;
        yo = y;
        zo = z;

    end

end

function J = jacobcalc(x, y, z)

    % Initialize Jacobian Matrix
    J = zeros(3,3); % Creates Jacobian that is 3 by 3 that has all zeros

    % Specify the delta values for Jacobian matrix
    dx = 1e-6 * x + eps;
    dy = 1e-6 * y + eps;
    dz = 1e-6 * z + eps;

    % Make function calls
    f = funcalc(x,y,z);
    fx = funcalc(x + dx,y,z);
    fy = funcalc(x,y + dy,z);
    fz = funcalc(x,y,z + dz);

    % Fill Jacobian Matrix
    % Row 1
    J(1,1) = (fx(1) - f(1))/dx;
    J(1,2) = (fy(1) - f(1))/dy;
    J(1,3) = (fz(1) - f(1))/dz;
    % Row 2
    J(2,1) = (fx(2) - f(2))/dx;
    J(2,2) = (fy(2) - f(2))/dy;
    J(2,3) = (fz(2) - f(2))/dz;
    % Row 3
    J(3,1) = (fx(3) - f(3))/dx;
    J(3,2) = (fy(3) - f(3))/dy; 
    J(3,2) = (fz(3) - f(3))/dz;
end

function f = funcalc(x,y,z)

    % Initialize function vector
    f = zeros(3,1);

    % Evaluate Function Vector

    f(1) = sin(x) + y.^2 + log(z) - 7;
    f(2) = 3*x + 2^y - z^3 + 1;
    f(3) = x + y + z -5;

end