function c = elasticityC2D(E,nu,TwoDType)
% ELASTICITYC2D Helper function to create plane-stress, plane-strain and axisymmetricc.
% Also works for function handle inputs
% Copyright 2022-2023 The MathWorks, Inc.

if ~isa(E,'function_handle') && ~isa(nu,'function_handle') % Retain original workflow for numeric E and nu
    G = E/(2*(1+nu));
    
    if strcmp(TwoDType,'ps')
        mu = 2*G*nu/(1-nu);
    else
        mu = 2*G*nu/(1-2*nu); % Both pn and as
    end
    
    c11 = [2*G+mu;0;0;G];
    c12 = [0;mu;G;0];
    c21 = [0;G;mu;0];
    c22 = [G;0;0;2*G+mu];
    c = [c11;c12;c21;c22];
    
    if strcmp(TwoDType,'as')
        c = @(region,state) 2*pi*c*region.x;
    end
else % E or nu or both are function handles
    if ~isa(E,'function_handle') % Convert numeric E to function handle
        E = @(l,s) E*ones(1,length(l.x));
    end

    if ~isa(nu,'function_handle')
        nu = @(l,s) nu*ones(1,length(l.x));
    end
    % Define all relevant quantities as in the numeric workflow as a function handles 
    G = @(l,s) E(l,s)./(2*(1+nu(l,s))); 
    
    if strcmp(TwoDType,'ps')
        mu = @(l,s) 2*G(l,s).*nu(l,s)./(1-nu(l,s));
    else
        mu = @(l,s) 2*G(l,s).*nu(l,s)./(1-2*nu(l,s)); % Both pn and as
    end
    c = @(l,s) pde.internal.FESolverUtilities.elasticityFcnHandleC2D(l,s,G,mu);
    
    if strcmp(TwoDType,'as')
        c = @(l,s) 2*pi*c(l,s).*l.x;
    end
end
end