Optim-nonlinear equation in matlab code

Question

I updated the question to clarify it more. Here is a graph:

For the curve in the attached photo, I hope to draw the curve. I have its equation and it is after simplification will be like this one

% Eq-2
    (b*Y* cos(v) + c    -  k*X*sin(v))^2 + ...
sqrt(k*X*(cos(v) + 1.0) +  b*Y*sin(v))^2) - d = 0.0 

Where:

v = atan((2.0*Y)/X) + c 

and b, c, d and k are constants.

from the attached graph,

The curve is identified in two points:

p1 @ (x=0)
p2 @ (y=0)

I a new on coding so accept my apologize if my question is not clear.

Thanks

Answer 1

So, I was having fun with this (thanks for that)!

Different question, different answer.

The solution below first searches for the constants causing the X and Y intercepts you were looking for (p1 and p2). For those constants that best fit the problem, it makes a plot, taking into account numerical issues.

In fact, you don't need eq. 1, because that's true always for any curve -- it's just there to confuse you, and problematic to use.

So, here it is:

function C = solve_for_F()

    % Points of interest
    px = 6;
    py = 4.2;

    % Wrapper function; search for those constants 
    % causing the correct X,Y intercepts (at px, py)
    G = @(C) abs(F( 0, px, C)) + ... % X intercept at px
             abs(F(py,  0, C));      % Y intercept at py

    % Initial estimate, based on your original equation
    C0 = [5/33
          1247745517111813/562949953421312
          4243112111277797/4503599627370496
          5/66];

    % Minimize the error in G by optimizing those constants
    C = fminsearch(G, C0);

    % Plot the solutions
    plot_XY(px, py, C);

end

function plot_XY(xmax,ymax, C)

    % graininess of X
    N = 100;

    % Find solutions for all alphae
    Y       = zeros(1,N);
    X       = linspace(0, xmax, N);
    y0      = linspace(ymax, 0, N);
    options = optimset('Display', 'off',...,...
                       'TolX'   , 1e-10);

    % Solve the nonlinear equation for each X
    for ii = 1:numel(X)

        % Wrapper function for fzero()
        fcn1 = @(y)F(y, X(ii), C);

        % fzero() is probably the fastest and most intuitive 
        % solver for this problem
        [Y(ii),~,flag] = fzero(fcn1, y0(ii), options);

        % However, it uses an algorithm that easily diverges
        % when the function slope is large. For those cases, 
        % solve with fminsearch()
        if flag ~= 1

            % In this case, the minimum of the absolute value
            % is searched for (which should be zero)
            fcn2 = @(y) abs(fcn1(y));

            Y(ii) = fminsearch(fcn2, y0(ii), options);
        end

    end

    % Now plot the X,Y solutions
    plot(X, Y,...
        'linewidth', 2,...
        'color',     [1 0.65 0]);    
    xlabel('X'), ylabel('Y')
    axis([0 xmax+.1 0 ymax+.1])

end

function fval = F(Y, X, C)

    % Unpack constants
    b = C(1);  d = C(3);
    c = C(2);  k = C(4);

    % pre-work
    V = atan2(2*Y, X) + c;

    % Eq. 2
    fval = sqrt( (b*Y*sin(V) + k*X*(cos(V) + 1))^2 + ...
                 (b*Y*cos(V) - k*X* sin(V)     )^2  ) - d;
end

Answer 2

So, after your edit, it is a bit more clear what you want.

I insist that your equation needs work -- the original equation (before your edit) simplified to what I have below. The curve for that looks like your plot, except the X and Y intercepts are at different locations, and funky stuff happens near X = 0 because you have numerical problems with the tangent (you might want to reformulate the problem).

But, after checking your equation, the following code should be helpful:

function solve_for_F()

    % graininess of alpha
    N = 100;

    % Find solutions for all alphae
    X       = zeros(1,N);
    options = optimset('Display', 'off');    
    alpha   = linspace(0, pi/2, N);
    x0      = linspace(6, 0, N);

    for ii = 1:numel(alpha)    
        X(ii) = fzero(@(x)F(x, alpha(ii)), x0(ii), options);
    end

    % Convert and make an X-Y plot
    Y = X .* tan(alpha);

    plot(X, Y,...
        'linewidth', 2,...
        'color',     [1 0.65 0]);

end

function fval = F(X, alpha)

    Y = X*tan(alpha);

    % Please, SIMPLIFY in the future
    A = 1247745517111813/562949953421312;
    B = 4243112111277797/4503599627370496;
    V = atan2(2*Y,X) + A;

    eq2 = sqrt(  (5/33*( Y*sin(V) + X/2*(cos(V) + 1) ))^2 + ...
                 (5/33*( Y*cos(V) - X/2* sin(V)      ))^2  ) - B;

    fval = eq2;

end

Results: