classdef (Sealed) StructuralMaterialAssignment < handle & matlab.mixin.CustomDisplay
% StructuralMaterialAssignment - Structural material property assignment over domain or subdomain
%  Structural analysis workflow in PDE Toolbox supports small deformation
%  linear elasticity with isotropic material. The properties that fully
%  characterize isotropic material are Young's modulus and Poisson's ratio.
%  In addition, mass density of the material is required to model body
%  loads, such as gravity. The StructuralMaterialAssignment object stores
%  the specified material properties on a domain or subdomain and appends
%  the object to the MaterialProperties property of the StructuralModel. In
%  multidomain geometric model, material properties can be assigned
%  independently for each subdomain or it can be one global assignment.
%
% Instances of this class can only be created by calling the
% structuralProperties method of the StructuralModel class.
%
% See also pde.StructuralModel, pde.StructuralModel/structuralProperties

% Copyright 2017-2021 The MathWorks, Inc.


properties (SetAccess = private)
    %RegionType - Type of geometric region to which the material properties are assigned
    % The type of geometric region specified as 'Cell' or 'Face'.
    RegionType;
end

properties
    % RegionID - IDs of the geometric regions to which the material properties are assigned
    % RegionID cannot exceed the number of faces for a 2-D geometry and the
    % number of cells for a 3-D geometry.
    RegionID;

    %YoungsModulus - Young's modulus of the material
    % A numeric scalar representing the Young's modulus or elastic
    % modulus of an isotropic material.
    YoungsModulus;

    %PoissonsRatio - Poisson's ratio of the material
    % A numeric scalar representing the Poisson's ratio of an isotropic
    % material.
    PoissonsRatio;
    
    %MassDensity  - Mass density of the material
    % A numeric scalar representing the mass density of an isotropic
    % material.
    MassDensity;
    
    %CTE - Coefficient of thermal expansion of the material
    % A numeric scalar representing the linear coefficient of thermal
    % expansion of an isotropic material.
    CTE;

    %HystereticDamping - Hysteretic (structural) damping coefficient
    HystereticDamping;
end

methods
    function self = StructuralMaterialAssignment(smar,varargin)
        self.RecordOwner = smar;
        parser = inputParser;
        parser.addParameter('face', []);
        parser.addParameter('cell', []);
        parser.addParameter('YoungsModulus', []);
        parser.addParameter('PoissonsRatio', []);
        parser.addParameter('MassDensity',[]);
        parser.addParameter('CTE',[]);
        parser.addParameter('HystereticDamping',[]);
        parser.parse(varargin{:});

        numdims = 2;
        if ~isempty(parser.Results.face)
            self.RegionType = 'Face';
            self.RegionID = parser.Results.face;
            systemsize = 2;
        elseif ~isempty(parser.Results.cell)
            self.RegionType = 'Cell';
            self.RegionID = parser.Results.cell;
            numdims = 3;
            systemsize = 3;
        end
        
        self.YoungsModulus = parser.Results.YoungsModulus;
        self.PoissonsRatio = parser.Results.PoissonsRatio;
        
        if ~isempty(parser.Results.MassDensity) % Mass density can be empty
            self.MassDensity = parser.Results.MassDensity;
        end
        
        if ~isempty(parser.Results.CTE) % CTE can be empty
            self.CTE = parser.Results.CTE;
        end
        
        if ~isempty(parser.Results.HystereticDamping) % HystereticDamping can be empty
            self.HystereticDamping = parser.Results.HystereticDamping;
        end

        self.checkAllMatrixMtlSizes(systemsize, numdims);
    end
end


methods % Setter methods

    function set.RegionType(self, rtype)
        self.RegionType = rtype;
    end
    
    function set.RegionID(self, rids)
        if ~isempty(self.RecordOwner) %#ok<MCSUP>
            nRegion =  self.RecordOwner.nRegion; %#ok<MCSUP>
            self.ValidateRegionID(rids,nRegion);
        end
        self.RegionID = rids;
    end
    
    function set.MassDensity(self, coef)
        self.CoefPrecheck(coef,'MassDensity');
        self.MassDensity = coef;
    end
    
    function set.YoungsModulus(self, coef)
        self.CoefPrecheck(coef,'YoungsModulus');
        self.YoungsModulus = coef;
    end
    
    function set.PoissonsRatio(self, coef)
        self.CoefPrecheck(coef,'PoissonsRatio');
        self.PoissonsRatio = coef;
    end

    function set.CTE(self, coef)
        %Negative and zero CTE are allowed.
        validateattributes(coef,{'double'},{'nonnan','nonempty','nonsparse','real','finite'},'CTE');
        self.CTE = coef;
    end

    function set.HystereticDamping(self, coef)
        validateattributes(coef,{'double'},{'nonnegative','nonnan','nonsparse','real','finite'},'HystereticDamping');
        self.HystereticDamping = coef;
    end
end


methods(Hidden=true, Access = {?pde.StructuralMaterialAssignmentRecords})
    function tf=sameCoefficients(self,other,varargin)
        if ~coefficientsMatch(self, other)
            tf = false;
            return
        end
        
        if isempty(varargin)
            tf = (isequal(self.MassDensity, other.MassDensity) && ...
                  isequal(self.YoungsModulus, other.YoungsModulus) && ...
                  isequal(self.PoissonsRatio, other.PoissonsRatio)&& ...
                  isequal(self.CTE, other.CTE)) && ...
                  isequal(self.HystereticDamping,other.HystereticDamping);
        else
            cname = varargin{1};
            switch cname
                case 'YoungsModulus'
                    tf = isequal(self.YoungsModulus, other.YoungsModulus);
                case 'PoissonsRatio'
                    tf = isequal(self.PoissonsRatio, other.PoissonsRatio);
                case 'MassDensity'
                    tf = isequal(self.MassDensity, other.MassDensity);
                case 'CTE'
                    tf = isequal(self.CTE, other.CTE);
                case 'HystereticDamping'
                    tf = isequal(self.HystereticDamping, other.HystereticDamping);
            end
        end
    end

    function tf=numericMaterialProperties(self,varargin)
        if isempty(varargin)
            tf = ~any(isa(self.YoungsModulus, 'function_handle')  || ...
                      isa(self.MassDensity,   'function_handle')  || ...
                      isa(self.PoissonsRatio, 'function_handle') || ...
                      isa(self.CTE, 'function_handle') || ...
                      isa(self.HystereticDamping, 'function_handle') );
        else
            cname = varargin{1};
            switch cname
                case 'YoungsModulus'
                    tf = ~isa(self.YoungsModulus, 'function_handle');
                case 'MassDensity'
                    tf = ~isa(self.MassDensity, 'function_handle');
                case 'PoissonsRatio'
                    tf = ~isa(self.PoissonsRatio, 'function_handle');
                case 'CTE'
                    tf = ~isa(self.CTE, 'function_handle');
                case 'HystereticDamping'
                    tf = ~isa(self.HystereticDamping, 'function_handle');
            end
        end
        
        %function_handle not supported for isotropic material.
        if ~tf
            error(message('pde:materialPropertySpecification:invalidMtlPropNonNumeric'));
        end        
        
    end


    function tf = coefficientsMatch(self, other)
        tf = true;
        tf = tf & (self.YoungsModulusDefined() == other.YoungsModulusDefined());
        tf = tf & (self.MassDensityDefined() == other.MassDensityDefined());
        tf = tf & (self.PoissonsRatioDefined() == other.PoissonsRatioDefined());
        tf = tf & (self.CTEDefined() == other.CTEDefined());
    end

    function performSolverPrecheck(self, systemsize, numfaces, numcells)
        ndims = 2;
        if strcmp(self.RegionType, 'Face')
            if any(self.RegionID > numfaces)
                error(message('pde:materialPropertySpecification:invalidFaceIndexPresolve'));
            end
        else
            if any(self.RegionID > numcells)
                error(message('pde:materialPropertySpecification:invalidCellIndexPresolve'));
            end
            ndims = 3;
        end
        checkAllMatrixMtlSizes(self, systemsize, ndims);
    end

end


methods(Hidden=true, Access = private)
    function checkAllMatrixMtlSizes(self, systemsize, ndims)
        if ~isempty(self.MassDensity)
            self.checkMassDensitySize(self.MassDensity, systemsize, ndims);
        end
        
        if ~isempty(self.CTE)
            self.checkMassDensitySize(self.CTE, systemsize, ndims);
        end
        
        if ~isempty(self.HystereticDamping)
            self.checkHystereticDampingSize(self.HystereticDamping, systemsize, ndims);
        end

        self.checkYoungsModulusSize(self.YoungsModulus, systemsize, ndims);
        self.checkPoissonsRatioSize(self.PoissonsRatio, systemsize, ndims);
    end
    %
end % Private methods



methods(Hidden=true, Access = {?pde.StructuralMaterialAssignmentRecords})
    function tf = YoungsModulusDefined(self)
        tf = self.coefDefined(self.YoungsModulus);
    end

    function tf = MassDensityDefined(self)
        tf = self.coefDefined(self.MassDensity);
    end

    function tf = PoissonsRatioDefined(self)
        tf = self.coefDefined(self.PoissonsRatio);
    end
    
    function tf = CTEDefined(self)
        tf = self.coefDefined(self.CTE);
    end
        
     function tf = HystereticDampingDefined(self)
        tf = self.coefDefined(self.HystereticDamping);
    end

    function tf = hasComplexCoefficient(self, loc, state)
        tf = false(4,1);
        tf(1) = pde.StructuralMaterialAssignment.coefIsComplexNumericOrFcnHdl(self.YoungsModulus, loc, state);
        tf(2) = pde.StructuralMaterialAssignment.coefIsComplexNumericOrFcnHdl(self.MassDensity, loc, state);
        tf(3) = pde.StructuralMaterialAssignment.coefIsComplexNumericOrFcnHdl(self.PoissonsRatio, loc, state);
        tf(4) = pde.StructuralMaterialAssignment.coefIsComplexNumericOrFcnHdl(self.CTE, loc, state);
        tf(5) = pde.StructuralMaterialAssignment.coefIsComplexNumericOrFcnHdl(self.HystereticDamping, loc, state);
        
        if any(tf)
            error(message('pde:materialPropertySpecification:invalidMtlPropValueComplex'));
        end
        
    end

end % end of methods


methods(Static, Access = private)
    function tf = coefDefined(coef)
        tf = (isnumeric(coef) && ~(isscalar(coef) && coef == 0) );
    end

    function tf = coefIsComplexNumericOrFcnHdl(coef, ~, ~)
        tf = false;
        if ~pde.StructuralMaterialAssignment.coefDefined(coef)
            return
        end
        if isnumeric(coef)
            if ~isreal(coef)
                tf = true;
            end
        end
    end

    
    function ok=ValidateRegionID(rgnid,nregion)
        % Must be real(non-complex), full, natural number.
        if ~isreal(rgnid) || ~all(rgnid(:) > 0) || issparse(rgnid) || any(mod(rgnid(:),1)~=0) || any(rgnid > nregion)
            error(message('pde:materialPropertySpecification:invalidRegionID'));
        end
        ok = true;
    end

    function ok=CoefPrecheck(coef,coefName)
        validateattributes(coef,{'double','function_handle'},{},coefName)
        if isa(coef,'function_handle')
            error(message('pde:materialPropertySpecification:invalidMtlPropNonNumeric'));
        else
            validateattributes(coef,{'double'},{'nonnan','nonempty','nonsparse','real','finite','nonnegative'},coefName);
            if all(coef==0) % Matrix form of coef can have some zeros, but not all.
                error(message('pde:materialPropertySpecification:invalidMtlPropNonPositive'));
            end
            
        end
        ok = true;
    end

    function checkMassDensitySize(coef, ~,~)
        if ~isscalar(coef)
            error(message('pde:materialPropertySpecification:invalidMtlPropNonNumeric'));
        end
    end

    function checkYoungsModulusSize(coef,~,~)
        if ~isscalar(coef)
            error(message('pde:materialPropertySpecification:invalidMtlPropNonNumeric'));
        end
    end
    
    function checkPoissonsRatioSize(coef,~,~)
        if ~isscalar(coef)
            error(message('pde:materialPropertySpecification:invalidMtlPropNonNumeric'));
        end
    end

    function checkCTESize(coef,~,~)
        if ~isscalar(coef)
            error(message('pde:materialPropertySpecification:invalidMtlPropNonNumeric'));
        end
    end

    function checkHystereticDampingSize(coef,~,~)
        if ~isscalar(coef)
            error(message('pde:materialPropertySpecification:invalidMtlPropNonNumeric'));
        end
    end

end

methods(Hidden=true)
    function delete(self)
        if ~isempty(self.RecordOwner) && isvalid(self.RecordOwner)
            self.RecordOwner.delistMaterialAssignments(self);
        end
    end
end


properties (Hidden = true, Access='private')
    RecordOwner;
end


end
