classdef (Sealed) ElectromagneticSourceAssignment  < pde.PDECompactDisplay
    % ElectromagneticSourceAssignment - Electromagnetic source assignment over domain or subdomain 
    % The ElectromagneticSourceAssignment object stores the specified
    % sources on a domain or subdomain and appends the object to the
    % Sources property of ElectromagneticModel. In a multidomain
    % geometric model, sources 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
    % electromagneticSource method of the ElectromagneticModel class.
    %
    % See also pde.ElectromagneticModel, pde.ElectromagneticModel/electromagneticSource
    
    % Copyright 2020-2022 The MathWorks, Inc.
    
    properties (SetAccess = private)
        %RegionType - Type of geometric region to which the electromagnetic
        %source is assigned.
        RegionType;
    end
    
    properties
        % RegionID - IDs of the geometric regions to which the electromagnetic
        % source is assigned.
        RegionID;
        
        %ChargeDensity - Electric charge density as source.
        ChargeDensity;
        
        %CurrentDensity - Electric current density as source.
        CurrentDensity;

        %Magnetization - Magnetization as a source.
        Magnetization;
    end
    
    methods
        function self = ElectromagneticSourceAssignment(emar,sourceInput)
            arguments
                emar pde.ElectromagneticSourceAssignmentRecords
                sourceInput struct
            end
            
            self.RecordOwner = emar;
            validDomain = emar.ParentModel.validDomain;
            self.RegionType = validDomain;
            self.RegionID = sourceInput.(validDomain);
            
            sourceInput = rmfield(sourceInput,validDomain);
            inputFields = fields(sourceInput);
            for f = 1:numel(inputFields)
                self.(inputFields{f}) = sourceInput.(inputFields{f});
            end
            self.PropertiesToValidate = setdiff(properties(self),{'RegionType','RegionID'});
        end
    end
    
    methods % Setter methods
        function set.RegionID(self, rids)
            ValidateRegionID(self,rids);
            self.RegionID = rids;
        end
        
        function set.ChargeDensity(self, coef)
            self.CoefPrecheck(coef,'ChargeDensity');
            self.ChargeDensity = coef;
        end
        
        function set.CurrentDensity(self, coef)
            self.CoefPrecheck(coef,'CurrentDensity');
            self.CurrentDensity = coef;
        end

        function set.Magnetization(self,coef)
            self.CoefPrecheck(coef,'Magnetization');
            self.Magnetization = coef;
        end
    end % end of setter methods
    
    methods(Hidden=true, Access = {?pde.ElectromagneticSourceAssignmentRecords})
        function tf=numericSource(self,varargin)
            tf = true;
            for prop = self.PropertiesToValidate.'
                tf = tf && ~isa(self.(prop{1}),'function_handle');
            end
        end
        
        function validateTypeSizeComplexity(self, loc, state, complexAllowed)
            for prop = self.PropertiesToValidate.'
                propName = prop{1};
                propValue = self.(prop{1});
                if isempty(propValue)
                    continue
                end

                if isa(propValue,'pde.ConductionResults')
                    validateConductionResults(self,propValue);
                    continue
                end
                
                coeffSize = allowedSize(self,propName);
                [complexityInput,sizeInput] = self.RecordOwner.ParentModel.validateFcnHdl(propValue, loc, state);

                if complexityInput && ~complexAllowed
                    regionInfo = pde.ElectromagneticModel.getRegionText(self.RegionType,self.RegionID);
                    error(message('pde:ElectromagneticModel:invalidComplex',propName,regionInfo));
                end

                if (sizeInput ~= coeffSize) && ~contains(self.RecordOwner.ParentModel.AnalysisType,'harmonic')
                    error(message('pde:ElectromagneticModel:incorrectSizeCoefFunction', propName, coeffSize));
                end
            end
        end
        
        function ok=ValidateRegionID(self,rgnid)
            if isempty(self.RecordOwner)
                return; % Skips validation of rgnid while loading model from file.
            end
            % Must be vector of natural numbers.
            maxDomainID = self.RecordOwner.ParentModel.maxDomainID;
            if ~isreal(rgnid) || ~all(rgnid(:) > 0) || issparse(rgnid) || any(mod(rgnid(:),1)~=0) || any(rgnid > maxDomainID)
                error(message('pde:materialPropertySpecification:invalidRegionID'));
            end
            ok = true;
        end
        
       function ok=CoefPrecheck(self,coef,coefName)
            if isempty(coef)
                return;
            end
            if isempty(self.RecordOwner)
                return; % Skips validation of rgnid while loading model from file.
            end
            
            % Validate that input is either numeric or function handle, or conduction results.
            validtypes = {'double','function_handle'};
            if self.RecordOwner.ParentModel.AnalysisType == "magnetostatic" && ~self.RecordOwner.ParentModel.IsTwoD
                % Validate that input can include pde.ConductionResults for
                % 3-D magnetostatic analysis
                validtypes = {'double','function_handle','pde.ConductionResults'};
            end
            validateattributes(coef,validtypes,{},'',coefName);

            if ~isnumeric(coef) 
                % Function handle and results objects are validated separately.
                ok = true;
                return
            end
            
            dim = allowedSize(self,coefName);

            if dim == 1
                shapeAttribute = 'scalar';
            else
                shapeAttribute = {'vector','numel',dim};
            end
            
            % validate const source based on type and dim
            if isequal(string(self.RecordOwner.ParentModel.AnalysisType),"harmonic")
                % the shape attribute has to be pushed to after solve to
                % know the appropriate size
                validateattributes(coef(:),{'double'},{'nonnan','nonempty','nonsparse','finite'},coefName);
            else
                validateattributes(coef(:),{'double'},[{'nonnan','nonempty','nonsparse','real','finite'} shapeAttribute(:)'],coefName);
            end
            ok = true;
       end
        
       function dim = allowedSize(self,coefName)
               % determine dimension for problem type
            dim = 1;
            if self.RecordOwner.ParentModel.IsTwoD && contains(self.RecordOwner.ParentModel.AnalysisType,'harmonic') ...
                    && isequal(string(self.RecordOwner.ParentModel.FieldType),"electric")
                dim = 2;
            elseif ~self.RecordOwner.ParentModel.IsTwoD ...
                    && ( isequal(coefName,"CurrentDensity") || isequal(string(self.RecordOwner.ParentModel.AnalysisType),"harmonic") )
                dim = 3;
            end
            if isequal(coefName,"Magnetization")
                if self.RecordOwner.ParentModel.IsTwoD
                    dim = 2;
                else
                    dim = 3;
                end
            end
       end


       function validateConductionResults(self,condRes)
           if isempty(condRes.CurrentDensity)
               error(message('pde:ElectromagneticModel:emptyConductionResults'));
           end

           %Check location and volume of conduction problem domain using mesh data.
           meshMS = self.RecordOwner.ParentModel.Mesh;
           meshConduction = condRes.Mesh;

           if ~isempty(meshMS)
               %First check bounding box.
               boundingBoxMS = [min(meshMS.Nodes,[],2); max(meshMS.Nodes,[],2)];
               boundingBoxConduction = [min(meshConduction.Nodes,[],2); max(meshConduction.Nodes,[],2)];
               if any(abs(boundingBoxMS - boundingBoxConduction) > eps(meshMS.MaxElementSize))
                   error(message('pde:ElectromagneticModel:mismatchedConductionResults'));
               end

               %Check volume
               volumeMS = volume(meshMS);
               volumeConduction = volume(meshConduction);
               if abs(volumeMS-volumeConduction) > eps(volumeMS)
                   error(message('pde:ElectromagneticModel:mismatchedConductionResults'));
               end
           end

       end
    end % end of methods that record keeper calls
    
    methods(Hidden=true)
        function delete(self)
            if ~isempty(self.RecordOwner) && isvalid(self.RecordOwner)
                self.RecordOwner.delistSourceAssignments(self);
            end
        end
    end


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