function ca = findStructuralIC(self, regiontype, regionid)  
%findStructuralIC - Find structural initial conditions assigned to a geometric region
% SICA = findStructuralIC(SICR, REGIONTYPE, REGIONID) returns a
% GeometricStructuralICs or NodalStructuralICs object, SICA, that defines
% the structural initial conditions specified for REGIONTYPE and REGIONID.
% Here, SICR is the InitialConditions property of the StructuralModel
% object. REGIONTYPE is 'Cell', 'Face', 'Edge', or 'Vertex'. REGIONID is
% the ID of the geometric region.
%  
% % Example: Find Initial Conditions for Cells of 3-D Geometry 
%
% % Create a geometry that consists of three nested cylinders, and include
% % the geometry in a structural dynamics model.
% gm = multicylinder([5 10 15],2);
% structuralmodel = createpde('structural','transient-solid');
% structuralmodel.Geometry = gm;
% pdegplot(structuralmodel,'CellLabels','on','FaceAlpha',0.5)
%
% % Set initial conditions for each cell.
% structuralIC(structuralmodel,'Displacement',[0,0,0],'Velocity',[0,0,0],'Cell',1);
% structuralIC(structuralmodel,'Displacement',[0,0.1,0],'Cell',2);
% structuralIC(structuralmodel,'Velocity',[0,0.2,0],'Cell',3);
%
% % Check the initial condition specification for cell 1.
% SICACell1 = findStructuralIC(structuralmodel.InitialConditions,'Cell',1)
%
% See also pde.StructuralModel/structuralIC

% Copyright 2017 The MathWorks, Inc.

narginchk(3,3)

if isempty(self.ParentPdemodel.Geometry)
    ca = [];
    return;
end

parser = inputParser;
parser.addParameter('Cell', [], @pde.EquationModel.isValidEntityID); 
parser.addParameter('Face', [], @pde.EquationModel.isValidEntityID);  
parser.addParameter('Edge', [], @pde.EquationModel.isValidEntityID);
parser.addParameter('Vertex', [], @pde.EquationModel.isValidEntityID);
parser.parse(regiontype, regionid);

if ~ismember('Vertex', parser.UsingDefaults)
   qregiontype = 'Vertex';
   qregionid = parser.Results.Vertex;
   maxid = self.ParentPdemodel.Geometry.NumVertices;   
elseif ~ismember('Edge', parser.UsingDefaults)
   qregiontype = 'Edge';
   qregionid = parser.Results.Edge;
   maxid = self.ParentPdemodel.Geometry.NumEdges;   
elseif ~ismember('Face', parser.UsingDefaults)
   qregiontype = 'Face';
   qregionid = parser.Results.Face;
   maxid = self.ParentPdemodel.Geometry.NumFaces;
else       
    qregiontype = 'Cell';
    qregionid = parser.Results.Cell;
    if self.ParentPdemodel.IsTwoD
        maxid = 0;
    else
        maxid = self.ParentPdemodel.Geometry.NumCells;
    end       
end
if any(qregionid > maxid)
       error(message('pde:pdeInitialConditions:invalidRegionID'));
end    
    

numqueries = numel(qregionid);
numassigns = numel(self.StructuralICAssignments);
ca = repmat(self.StructuralICAssignments(1),size(qregionid));
for i = 1:numqueries
    rid = qregionid(i);
    icfound = false;
    for j = numassigns:-1:1
       thisic =  self.StructuralICAssignments(j);
        if isa(thisic, 'pde.NodalStructuralICs') || (strcmp(qregiontype,thisic.RegionType) && ismember(rid, thisic.RegionID))        
          ca(i) = thisic;
          icfound = true;
          break
       end
    end
    if ~icfound
      error(message('pde:structuralInitialConditions:entitiesWithoutICs'));
    end   
end

end