function [xC,Resnorm,FVAL,EXITFLAG,OUTPUT,LAMBDA,JACOB] = ...
    lsqncommon(funfcn,xC,lb,ub,options,defaultopt,optimgetFlag,caller,initVals, ...
    sizes,flags,mtxmpy,varargin)
%

%LSQNCOMMON Solves non-linear least squares problems.
%   It contains all the setup code common to both LSQNONLIN and LSQCURVEFIT to 
%   call the different algorithms.

%   Copyright 1990-2022 The MathWorks, Inc.

[sizes.fRows,sizes.fCols] = size(initVals.F);
sizes.nFun = numel(initVals.F);
F = initVals.F(:);
Jac = initVals.J;

hasJacobMult = false;
Jinfo = [];
if flags.grad
    % check Jacobian mutliply function
    hasJacobMult = ~isempty(optimget(options,'JacobMult',defaultopt,optimgetFlag));
    if hasJacobMult
        % User has a JacobianMultiplyFcn - the initial Jacobian is actually
        % "Jinfo", which can be anything. In order to bypass validation for
        % Jac, swap the "pointers".
        Jinfo = Jac;
        Jac = [];
    end
end


% Check to see if the function vector or the user-supplied Jacobian at the initial point are well-defined.
% If not, then terminate immediately.
if ~allfinite(F) 
    error('optimlib:commonMsgs:UsrObjUndefAtX0', ...
        getString(message('optimlib:commonMsgs:UsrObjInfOrNaNAtX0', caller)));    
end
if ~allfinite(Jac)
    error('optimlib:commonMsgs:UserJacUndefAtX0', ...
        getString(message('optimlib:commonMsgs:JacInfOrNaNAtX0',caller)));
end

% Check for non-double data typed values returned by user functions
if ~isempty( isoptimargdbl(caller, {'F','J'}, F, Jac) )
    error('optimlib:commonMsgs:NonDoubleFunVal',getString(message('optimlib:commonMsgs:NonDoubleFunVal',caller)));
end

% Populate finDiffFlags hasLBs and hasUBs - these must be set regardless of
% whether finite-differences are used since the algorithms rely on them
% being there.
% NOTE: we assume lb and ub are column vectors. checkbounds() performs this
% action and must be called before reaching here
finDiffFlags.hasLBs = isfinite(lb);
finDiffFlags.hasUBs = isfinite(ub);
finDiffFlags.anyFiniteBounds = any(finDiffFlags.hasLBs | finDiffFlags.hasUBs);

if finDiffFlags.anyFiniteBounds
    % Check for complex x0 and finite bounds
    if ~isreal(xC)
        error('optimlib:lsqncommon:BoundsAndComplexX0',getString(message('optimlib:lsqncommon:BoundsAndComplexX0')));
    end

    % Check for complex f(x0) or f'(x0) and finite bounds
    if ~isreal(F) || ~isreal(Jac)
        error('optimlib:lsqncommon:BoundsAndComplexInitialFunctionEval',...
        getString(message('optimlib:lsqncommon:BoundsAndComplexInitialFunctionEval')));
    end
end

if flags.grad
    % check size of Jacobian
    [Jrows, Jcols] = size(Jac);
    if ~hasJacobMult 
        % Not using 'JacobMult' so Jacobian must be correct size
        if Jrows ~= sizes.nFun || Jcols ~= sizes.nVar
            error(message('optimlib:lsqncommon:InvalidJacSize', sizes.nFun, sizes.nVar))
        end
    end
end

algorithm = optimget(options,'Algorithm',defaultopt,optimgetFlag);
initDamping = optimget(options,'InitDamping',defaultopt,optimgetFlag);
if ~iscell(algorithm)    
    initLMparam = initDamping; 
else
    initLMparam = algorithm{2}; % Initial Levenberg-Marquardt parameter
    algorithm = algorithm{1};   % Algorithm string
end

switch algorithm
    case 'trust-region-reflective'
        % trust-region-reflective and not enough equations -- switch to Levenberg-Marquardt
        if sizes.nFun < sizes.nVar
            warning(message('optimlib:lsqncommon:SwitchToLineSearch'))
            trustRegion = false; % Call levenbergMarquardt
        else
            trustRegion = true;
        end
    case 'levenberg-marquardt'
        trustRegion = false;
    otherwise
        error(message('optimlib:lsqncommon:InvalidAlgorithm', upper( caller )))
end


if strcmpi(optimget(options,'Diagnostics',defaultopt,optimgetFlag),'on')
    % Do diagnostics on information so far
    if trustRegion
        OUTPUT.algorithm = 'trust-region-reflective';
    else
        OUTPUT.algorithm = 'levenberg-marquardt';
    end
    diagnose(caller,OUTPUT,flags.grad,false,false,false,...
        xC(:),0,0,0,0,lb,ub,funfcn,{''});
end

% Prepare options, sizes and flags for finitedifferences
 % No nonlinear constraints
sizes.xShape = size(xC);
sizes.mNonlinEq = 0; 
sizes.mNonlinIneq = 0;

% Read options for finitedifferences - note some flags (hasLBs and hasUBs -
% defined above) must be set regardless since the algorithms rely on them
% being there.
useParallel = optimget(options,'UseParallel',defaultopt,optimgetFlag);
options.UseParallel = validateopts_UseParallel(useParallel,true,true);
DerivativeCheck = strcmp(optimget(options,'DerivativeCheck',defaultopt,optimgetFlag),'on');
if ~flags.grad || DerivativeCheck
    % Read in and error check option TypicalX
    [typicalx,ME] = getNumericOrStringFieldValue('TypicalX','ones(numberOfVariables,1)', ...
        ones(sizes.nVar,1),'a numeric value',options,defaultopt);
    if ~isempty(ME)
        throw(ME)
    end
    checkoptionsize('TypicalX', size(typicalx), sizes.nVar);
    options.TypicalX = typicalx;
    options.DiffMinChange = optimget(options,'DiffMinChange',defaultopt,optimgetFlag);
    options.DiffMaxChange = optimget(options,'DiffMaxChange',defaultopt,optimgetFlag);
    options.FinDiffType = optimget(options,'FinDiffType',defaultopt,optimgetFlag);
    options = validateFinDiffRelStep(sizes.nVar,options,defaultopt,optimgetFlag);

    % Create structure of flags for finitedifferences
    finDiffFlags.fwdFinDiff = strcmpi(options.FinDiffType,'forward'); % Check for forward fin-diff
    finDiffFlags.scaleObjConstr = false; % No scaling
    finDiffFlags.chkComplexObj = false;  % Don't check whether objective function values are complex
    finDiffFlags.isGrad = false;         % Compute Jacobian, not gradient

    % finDiffFlags.chkFunEval will be first used by validateFirstDerivatives().
    if trustRegion
        % Don't check function values during finite differences
        finDiffFlags.chkFunEval = false;
    else
        % Check function values during finite differences
        finDiffFlags.chkFunEval = true;
    end
end

% Already checked options.GradObj in lsqnsetup - just overwrite GradObj
if flags.grad
    options.GradObj = 'on';
end
options.GradConstr = 'off'; 

if flags.grad && DerivativeCheck
% Check derivatives
    validateFirstDerivatives(funfcn,{''},xC, ...
        lb,ub,options,finDiffFlags,sizes,varargin{:});
end

% Execute algorithm
if trustRegion
    if ~flags.grad % provide sparsity of Jacobian if not provided.
        Jstr = optimget(options,'JacobPattern',defaultopt,optimgetFlag);
        if isempty(Jstr) || ischar(Jstr) 
            % options.JacobPattern is the default or empty: 'sparse(ones(jrows,jcols))'
            % Bypass column re-ordering, graph coloring, etc. since we will
            % perform dense finite-differencing anyway.
            Jstr = []; 
            finDiffFlags.sparseFinDiff = false;
        else 
            % Pattern matrix  - other datatypes (cell-array, struct) are checked in optimset and its
            % helper functions
            checkoptionsize('JacobPattern', size(Jstr), sizes.nVar, sizes.nFun);
            finDiffFlags.sparseFinDiff = true;
        end
    else
        Jstr = [];
        if hasJacobMult
            % Now that validation is complete, restore Jinfo into Jac since
            % snls* will use it as Jinfo for the JacobMult.
            Jac = Jinfo;
        end
    end
    % Set MaxFunEvals appropriately for trust-region-reflective
    defaultopt.MaxFunEvals = '100*numberOfVariables';
    
    if any(lb==ub)
        [xC,FVAL,LAMBDA,JACOB,EXITFLAG,OUTPUT,msgData]= ...
            snlsFixedVar(funfcn,xC,lb,ub,flags.verbosity,options,defaultopt,F,Jac,caller, ...
            Jstr,flags.computeLambda,mtxmpy,flags.detailedExitMsg,finDiffFlags,varargin{:});
    else
        [xC,FVAL,LAMBDA,JACOB,EXITFLAG,OUTPUT,msgData]= ...
            snls(funfcn,xC,lb,ub,flags.verbosity,options,defaultopt,F,Jac,caller, ...
            Jstr,flags.computeLambda,mtxmpy,flags.detailedExitMsg,finDiffFlags,varargin{:});
    end

else
    % Set MaxFunEvals appropriately for Levenberg-Marquardt
    defaultopt.MaxFunEvals = '200*numberOfVariables';
    
    [xC,FVAL,LAMBDA,JACOB,EXITFLAG,OUTPUT,msgData] = ...
        levenbergMarquardt(funfcn,xC,lb,ub,flags.verbosity,options,defaultopt,F,Jac, ...
        caller,initLMparam,flags.detailedExitMsg,optimgetFlag,finDiffFlags,sizes,varargin{:});
end

Resnorm = FVAL'*FVAL;

ProblemdefOptions = optimget(options, 'ProblemdefOptions',defaultopt,optimgetFlag);
FromEqnSolve = false;
if ~isempty(ProblemdefOptions) && isfield(ProblemdefOptions, 'FromEqnSolve')
    FromEqnSolve = ProblemdefOptions.FromEqnSolve;
end

% Handle exit message creation
if FromEqnSolve 
    % If called from EquationProblem/solve, store msgData and defer
    % creation of the exit message.
    OUTPUT.msgData = msgData;
elseif flags.makeExitMsg
    % Check whether we need the "link-free" version of the message - which
    % saves the most costly part of createExitMsg, which is stripping the
    % hyperlinks from the message.
    if flags.makeOutputStruct
        % Call createExitMsg to form and print exit message, if requested.
        OUTPUT.message =  createExitMsg(msgData{:});
    else
        % Call createExitMsg just to print exit message.
        createExitMsg(msgData{:});
    end
end

% Reset FVAL to original shapes
FVAL = reshape(FVAL,sizes.fRows,sizes.fCols);

