FMINCON в функции ограничения есть дополнительные параметры, которых нет в целевой функции

Здравствуйте, я хочу решить проблему оптимизации массы конструкции - арки - с ограничением прочности. Я новичок в Matlab. Я пытаюсь сделать все с вложенными функциями в одном файле, потому что моя цел...

Здравствуйте, я хочу решить проблему оптимизации массы конструкции - арки - с ограничением прочности. Я новичок в Matlab. Я пытаюсь сделать все с вложенными функциями в одном файле, потому что моя целевая функция имеет две переменные, а функция ограничения имеет четыре переменные. Но, конечно, не получается. Не думаю, что это единственная проблема кода.

Я прошу прощения за мою глупость и большое спасибо за любую помощь!

Весь мой код:

function [x, fval] = optim1(x0)

fun = @obj;
x0 = [0.8, 0.2, 3, 5];
A = [1,-8,0,0;-1,3,0,0];
b = [0;0];
Aeq = [];
beq = [];
lb = [0.5,0.1,1,0];
ub = [1.5,0.6,6,30];
nonlcon = @stress;      

[x,fval] = fmincon(fun,x0,A,b,Aeq,beq,lb,ub,nonlcon);
    function [M] = obj (x1,x2)
                %% VOLUME
                z = sym ('z');   
                %% INTEGRAL - area S1
                    f = 6.5;
                    L = 30;
                    f1 = f + 0.5.*x1;
                    L1 = L + 0.5.*x1;
                    x_1 = (-0.5).*x1;
                        y1 = f1.*sqrt(L1.^2 - 4.*z.^2)./L1;
                        S1 = int(y1, z, x_1, L1);
                %% INTEGRAL - area S2
                    f = 6.5;
                    L = 30;
                    f2 = f - 0.5.*x1;
                    L2 = L - 0.5.*x1;
                    x_2 = 0.5.*x1;
                        y2 = f2.*sqrt(L2.^2 - 4.*z.^2)./L2;
                        S2 = int(y2, z, x_2, L2);
                %% Volume
                    v = (S1 - S2).*x2;
                %% Massa
                po = 0.5; % kg/m3
            M = v*po;
     end

    function [c,ceq] = stress (x1,x2,x3,x4)
                %% Data
                z = sym ('z');
                f = 6.5;
                L = 30;
                q_sv = 0.927;
                Rc = 11000;

                yk = f.*sqrt(L^2 + 4.*x4.^2); % arch axis
                a = atan((-4)*f.*z./(L*sqrt(L^2 - 4.*z.^2))); % for INTEGR (z): angle alfa
                a_s = atan((-4)*f.*x4./(L*sqrt(L^2 - 4.*x4.^2)));% for internal forces (х4 - variable): angle alfa
                phi = pi - a_s; % angle phi
                F = x2 .* x1; % AREA of section
                W = x2 .* (x1.^2)./6; % moment of resistance for section

                q_s = 5.355*cos(1.5.*a).*x3; % SNOW load
                T = taylor(q_s,z); %(z for integr)
                q_pokr = 0.309.*x3; % COATING load 
                %% REACTIONS
                % snow load integrals
                q_s1 = int(T, z, -0.98*L/2, -0.00001) + int(T, z, 0.00001, 0.98*L/2);
                q_s2 = int(T, z, -0.99*L/2, -0.00001);
                q_s3 = int(T, z, -0.99*L/2, x4);

                % reactions
                Va = ((q_sv + q_pokr)*(L^2)/2 + q_s1.*L/2)./L;
                H = (Va.*L/2 - (q_sv + q_pokr)*(L^2)/8 - q_s2 .* L/3)./f;

                %% INTERNAL forces
                Mk = H.*yk - Va.*x4 + (q_sv + q_pokr).*(x4.^2)./2 + q_s3.*(x4.^2)/2;
                Nk = sin(phi).*(q_sv + q_pokr).*x4 + cos(phi).*q_s3.*x4 - H.*cos(phi)-Va.*sin(phi);

                Md = Mk./(1-(Nk.*((66.89.*x1).^2)/(3000.*x2.*x1.*Rc)));
            str = (Nk/F) + (Md/W) - Rc;
            
        %% NONLCON
        c = [];
        ceq = str(x1,x2,x3,x4);
    end

end