function pdetool_solve(flag) %#ok flag is used as a switch using nargin
% case: solve PDE

%   Copyright 1994-2023 The MathWorks, Inc.

    pde_fig=findobj(allchild(0),'flat','Tag','PDETool');

    % if no geometry, create a default L-shape:
    pdegd=get(findobj(get(pde_fig,'Children'),'flat',...
        'Tag','PDEMeshMenu'),'UserData');
    if isempty(pdegd), pdetool('membrane'), end

    zstate=findobj(get(pde_fig,'Children'),'Tag','zoomBtn');
    if (zstate.UserData==1)
        pdezoom(pde_fig,'off')
        zstate.UserData = 0;
        opthndl=findobj(get(pde_fig,'Children'),'flat','Tag','PDEOptMenu');
        set(findobj(get(opthndl,'Children'),'flat','Tag','PDEZoom'),...
            'Checked','off')
    end

    flg_hndl=findobj(get(pde_fig,'Children'),'flat','Tag','PDEFileMenu');
    flags=get(flg_hndl,'UserData');
    flag1=flags(3);
    if abs(flag1)
        pdetool('changemode',0)
    end
    flags=get(flg_hndl,'UserData');
    if flags(3)==-1
        % error in decsg
        return;
    elseif ~flags(2)
        pdetool('cleanup')
    end
    flags=get(flg_hndl,'UserData');
    flag2=flags(4);
    if flag2
        pdetool('initmesh')
    end
    flags=get(flg_hndl,'UserData');
    oldmode=flags(2);
    flags(2)=3;                           % mode_flag=3
    set(flg_hndl,'UserData',flags)

    pdeinfo('Solving PDE...');
    set(pde_fig,'Pointer','watch')

    drawnow

    bndhndl=findobj(get(pde_fig,'Children'),'flat','Tag','PDEBoundMenu');
    bl=get(findobj(get(bndhndl,'Children'),'flat',...
        'Tag','PDEBoundMode'),'UserData'); %#ok - used in eval
    dl1=getappdata(pde_fig,'dl1'); %#ok - used in eval

    % Unpack parameters:
    params=get(findobj(get(pde_fig,'Children'),'flat','Tag','PDEPDEMenu'),...
        'UserData');
    ns=getappdata(pde_fig,'ncafd');
    nc=ns(1); na=ns(2); nf=ns(3); nd=ns(4);
    c=params(1:nc,:); %#ok - used in eval
    a=params(nc+1:nc+na,:); %#ok - used in eval
    f=params(nc+na+1:nc+na+nf,:); %#ok - used in eval
    d=params(nc+na+nf+1:nc+na+nf+nd,:); %#ok - used in eval

    pde_type=get(findobj(get(pde_fig,'Children'),'flat','Tag','PDEHelpMenu'),...
        'UserData');

    if pde_type>1
        timepar=getappdata(pde_fig,'timeeigparam');
        tlist=str2num(deblank(timepar(1,:))); %#ok - str2num may be ok here
        u0=deblank(timepar(2,:)); %#ok - used in eval
        ut0=deblank(timepar(3,:)); %#ok - used in eval
        r=str2num(deblank(timepar(4,:)));  %#ok - str2num may be ok here
        rtol=str2num(deblank(timepar(5,:))); %#ok - str2num may be ok here
        atol=str2num(deblank(timepar(6,:))); %#ok - str2num may be ok here
    end

    h=findobj(get(pde_fig,'Children'),'flat','Tag','PDEMeshMenu');
    hp=findobj(get(h,'Children'),'flat','Tag','PDEInitMesh');
    p=get(hp,'UserData');
    he=findobj(get(h,'Children'),'flat','Tag','PDERefine');
    e=get(he,'UserData');
    ht=findobj(get(h,'Children'),'flat','Tag','PDEMeshParam');
    t=get(ht,'UserData');

    solveparams=getappdata(pde_fig,'solveparam');
    adaptflag=str2double(deblank(solveparams(1,:)));
    nonlinflag=str2double(deblank(solveparams(7,:)));
    nltol=str2num(deblank(solveparams(8,:))); %#ok - str2num may be ok here
    nlinit=deblank(solveparams(9,:));
    jac=deblank(solveparams(10,:)); %#ok - used in eval
    nlinnorm=lower(deblank(solveparams(11,:)));
    if ~strcmp(nlinnorm,'energy')
        nlinnorm=str2num(nlinnorm);  %#ok - str2num may be ok here
    end

    err = false;
    errstr = '';
    % Solve PDE and catch error:
    if pde_type==1
        % solve elliptic problem:
        if adaptflag
            maxt=str2num(deblank(solveparams(2,:))); %#ok - str2num may be ok here
            ngen=str2num(deblank(solveparams(3,:))); %#ok - str2num may be ok here
            tselmet=deblank(solveparams(4,:)); %#ok - used in eval
            par=str2num(deblank(solveparams(5,:))); %#ok - str2num may be ok here
            refmet=deblank(solveparams(6,:)); %#ok - used in eval

            if ~nonlinflag
                nlin='off'; %#ok - used in eval
            else
                nlin='on'; %#ok - used in eval
            end
            if isempty(nlinit)
                try
                    eval(['[u,p,e,t]=adaptmesh(dl1,bl,c,a,f,''mesh'',p,e,t,',...
                        '''par'',par,''ngen'',ngen,''maxt'',maxt,''nonlin'',nlin,',...
                        '''toln'',nltol,''tripick'',tselmet,''rmethod'',refmet,',...
                        '''jac'',jac,''norm'',nlinnorm);']);
                catch ME
                    err = true;
                    errstr = ME.message;
                end
            else
                try
                    eval(['[u,p,e,t]=adaptmesh(dl1,bl,c,a,f,''mesh'',p,e,t,',...
                        '''par'',par,''ngen'',ngen,''maxt'',maxt,''nonlin'',nlin,',...
                        '''toln'',nltol,''init'',nlinit,''tripick'',tselmet,',...
                        '''rmethod'',refmet,''jac'',jac,''norm'',nlinnorm);']);
                catch ME
                    err = true;
                    errstr = ME.message;
                end
            end
            if ~err
                set(hp,'UserData',p), set(he,'UserData',e)
                set(ht,'UserData',t)
                % Update max no of triangle default for adaptive solver
                % Double the current no, or 1000, whichever is greatest:
                tristr=int2str(max(1.5*size(t,2),1000));
                setappdata(pde_fig,'solveparam',...
                    char(solveparams(1,:),tristr,solveparams(3:11,:)))

                meshstat=getappdata(pde_fig,'meshstat');
                meshstat=[meshstat 4];
                n=length(meshstat);
                setappdata(pde_fig,'meshstat',meshstat);

                % save this generation of the mesh for unrefine purposes
                setappdata(pde_fig,['p' int2str(n)],p);
                setappdata(pde_fig,['e' int2str(n)],e);
                setappdata(pde_fig,['t' int2str(n)],t);

                % Enable 'undo mesh change' option:
                meshhndl=findobj(get(pde_fig,'Children'),'flat','Tag','PDEMeshMenu');
                h=findobj(get(meshhndl,'Children'),'flat','Tag','PDEUnrefine');
                set(h,'Enable','on');
            end
        elseif nonlinflag
            if isempty(nlinit)
                try
                    eval(['[u,res]=pdenonlin(bl,p,e,t,c,a,f,''jacobian'',jac,',...
                        '''tol'',nltol,''norm'',nlinnorm,''report'',''on'');']);
                catch ME
                    err = true;
                    errstr = ME.message;
                end

            else
                try
                    eval(['[u,res]=pdenonlin(bl,p,e,t,c,a,f,''jacobian'',jac,',...
                        '''tol'',nltol,''U0'',nlinit,''norm'',nlinnorm,',...
                        '''report'',''on'');']);
                catch ME
                    err = true;
                    errstr = ME.message;
                end
            end
        else
            try
                eval('u=assempde(bl,p,e,t,c,a,f);');
            catch ME
                err = true;
                errstr = ME.message;
            end
        end
        if ~err
            if any(isnan(u)) && nonlinflag
                % u contains NaN's if problem is nonlinear and assembled
                % using 'assempde'
                pdeinfo('Switching to nonlinear solver...',0);
                nltol=str2num(deblank(solveparams(8,:))); %#ok - str2num may be ok here
                try
                    eval('[u,res]=pdenonlin(bl,p,e,t,c,a,f,nltol);');
                catch ME
                    err = true;
                    errstr = ME.message;
                end
            end
            if any(isnan(u))
                % if u still contains NaN's, problem is probably time dependent
                pdetool('error','Unable to solve elliptic PDE. Problem may be time dependent or nonlinear.')
                set(pde_fig,'Pointer','arrow')
                drawnow
                pdeinfo;
                return;
            end
        end
        l=[];
        pdeinfo('PDE solution computed.');
    elseif pde_type==2
        % solve parabolic problem:
        try
            eval('u=parabolic(u0,tlist,bl,p,e,t,c,a,f,d,rtol,atol);');
        catch ME
            err = true;
            errstr = ME.message;
        end
        if(~err && (size(u, 2) ~= size(tlist,2)))
            warndlg(getString(message('pde:pdetool:incompSoln')), ...
                getString(message('pde:pdetool:incompSolnTitle')), 'non-modal');
        end
        l=[];
        pdeinfo('PDE solution computed.');
    elseif pde_type==3
        % solve hyperbolic problem:
        try
            eval('u=hyperbolic(u0,ut0,tlist,bl,p,e,t,c,a,f,d,rtol,atol);');
        catch ME
            err = true;
            errstr = ME.message;
        end
        if(~err && size(u, 2) ~= size(tlist,2))
            warndlg(getString(message('pde:pdetool:incompSoln')), ...
                getString(message('pde:pdetool:incompSolnTitle')), 'non-modal');
        end
        l=[];
        pdeinfo('PDE solution computed.');
    elseif pde_type==4
        % solve eigenvalue problem:
        try
            eval('[u,l]=pdeeig(bl,p,e,t,c,a,d,r);');
        catch ME
            err = true;
            errstr = ME.message;
        end

        if ~err
            if isempty(l), %#ok - from eval
                plthndl=findobj(get(pde_fig,'Children'),'flat','Tag','PDEPlotMenu');
                set(plthndl,'UserData',u);
                winhndl=findobj(get(pde_fig,'Children'),'flat','Tag','winmenu');
                set(winhndl,'UserData',l);
                pdetool('error',' No eigenvalues in specified range.')
                set(pde_fig,'Pointer','arrow')
                drawnow
                pdeinfo('PDE solution computed.');
            else
                pdeinfo(sprintf('%i eigenvalues found. Use Plot Selection dialog box to select higher eigenmodes.',length(l)));
            end
        end
    end

    % Check if error:
    if err
        % Remove disturbing newlines:
        nl=find(abs(errstr)==10);
        if ~isempty(nl)
            rmnl=[];
            if nl(1)==1
                rmnl=1;
            end
            nonl=length(nl);
            if nonl>1
                if nl(nonl)==length(errstr) && nl(nonl-1)==length(errstr)-1
                    rmnl=[rmnl nl(nonl)];
                end
            end
            errstr(rmnl)=[];
        end
        pdetool('error',errstr);
        set(pde_fig,'Pointer','arrow')
        drawnow
        % Restore old mode
        flags=get(flg_hndl,'UserData');
        flags(2)=oldmode;
        set(flg_hndl,'UserData',flags)
        return;
    end

    % Save solution:
    plothndl=findobj(get(pde_fig,'Children'),'flat','Tag','PDEPlotMenu');
    set(plothndl,'UserData',u);
    % save eigenvalues:
    winhndl=findobj(get(pde_fig,'Children'),'flat','Tag','winmenu');
    set(winhndl,'UserData',l);

    % Enable export:
    solvehndl=findobj(get(pde_fig,'Children'),'flat','Tag','PDESolveMenu');
    set(findobj(get(solvehndl,'Children'),'flat','Tag','PDEExpSol'),...
        'Enable','on')

    % Set flags
    flags=get(flg_hndl,'UserData');
    flags(5)=0; flags(6)=1; flags(7)=0;   % flag3=0, flag4=1, flag5=0
    set(flg_hndl,'UserData',flags)

    plotflags=getappdata(pde_fig,'plotflags');
    if pde_type==2 || pde_type==3
        plotflags(12)=size(u,2);
    else
        plotflags(12)=1;
    end
    setappdata(pde_fig,'plotflags',plotflags)

    % Turn off replay of movie
    animparams=getappdata(pde_fig,'animparam');
    animparams(3)=0;
    setappdata(pde_fig,'animparam',animparams)

    % Update plot dialog box:
    pdeptdlg('initialize',1,getappdata(pde_fig,'plotstrings'));

    % flag is set if we're solving from PDEPTDLG. The solution plot
    % will be handled from PDEPTDLG.
    if plotflags(8) && ~isempty(u) && nargin==0
        % do plot solution automatically:
        pdeptdlg('plot')
    elseif adaptflag && (oldmode==2) && pde_type==1
        % We're still displaying the old mesh; update it
        pdetool meshmode
    else
        % Restore old mode since we are not plotting solution nor updating:
        flags=get(flg_hndl,'UserData');
        flags(2)=oldmode;
        set(flg_hndl,'UserData',flags)
    end

    set(pde_fig,'Pointer','arrow')
    drawnow

    if pde_type~=4
        pdeinfo('Select a new plot, or change mode to alter PDE, mesh, or boundaries.');
    end
end
