Index: /issm/trunk-jpl/src/m/classes/qmu.m =================================================================== --- /issm/trunk-jpl/src/m/classes/qmu.m (revision 24172) +++ /issm/trunk-jpl/src/m/classes/qmu.m (revision 24173) @@ -66,6 +66,6 @@ end end - if ~isempty(md.qmu.vpartition), - if numel(md.qmu.vpartition)~=md.mesh.numberofvertices + if ~isempty(md.qmu.vpartition) & ~any(isnan(md.qmu.vpartition)), + if numel(md.qmu.vpartition)~=md.mesh.numberofvertices md = checkmessage(md,['user supplied vertex partition for qmu analysis should have size md.mesh.numberofvertices x 1']); end @@ -77,5 +77,5 @@ end end - if ~isempty(md.qmu.epartition), + if ~isempty(md.qmu.epartition) & ~any(isnan(md.qmu.epartition)), if numel(md.qmu.epartition)~=md.mesh.numberofelements, md = checkmessage(md,['user supplied element partition for qmu analysis should have size md.mesh.numberofelements x 1']); @@ -88,7 +88,8 @@ end end - if isempty(md.qmu.vpartition) & isempty(md.qmu.epartition), - md = checkmessage(md,['for qmu analysis, at least one of element and partitioning vectors need to be supplied!']); - end + if isempty(md.qmu.vpartition) | any(isnan(md.qmu.vpartition)) | isempty(md.qmu.epartition) | any(isnan(md.qmu.epartition)), + md = checkmessage(md,['for qmu analysis, both an element and partitioning vectors need to be supplied with no nan values! One can be defaulted to all zeros.']); + end + end % }}} function disp(self) % {{{ Index: /issm/trunk-jpl/src/m/classes/qmu.py =================================================================== --- /issm/trunk-jpl/src/m/classes/qmu.py (revision 24172) +++ /issm/trunk-jpl/src/m/classes/qmu.py (revision 24173) @@ -25,5 +25,6 @@ self.params = OrderedStruct() self.results = OrderedDict() - self.partition = float('NaN') + self.vpartition = float('NaN') + self.epartition = float('NaN') self.numberofpartitions = 0 self.numberofresponses = 0 @@ -107,5 +108,6 @@ s+=" %-*s: [%ix%i] '%s'\n" % (maxlen+1,fname,a,b,type(getattr(result,fname))) - s+="%s\n" % fielddisplay(self,'partition','user provided mesh partitioning, defaults to metis if not specified') + s+="%s\n" % fielddisplay(self,'vpartition','user provided mesh partitioning (vertex based)') + s+="%s\n" % fielddisplay(self,'epartition','user provided mesh partitioning (element based)') s+="%s\n" % fielddisplay(self,'numberofpartitions','number of partitions for semi-discrete qmu') s+="%s\n" % fielddisplay(self,'variabledescriptors','') @@ -152,11 +154,22 @@ md.checkmessage('in parallel library mode, the requirement is for md.cluster.np = md.qmu.params.processors_per_evaluation * number_of_slaves, where number_of_slaves will automatically be determined by Dakota. Modify md.cluster.np accordingly') - if np.size(md.qmu.partition) > 0: - if (np.size(md.qmu.partition)!=md.mesh.numberofvertices and np.size(md.qmu.partition) != md.mesh.numberofelements) or np.size(md.qmu.partition,1) != 1: - md.checkmessage("user supplied partition for qmu analysis should have size (md.mesh.numberofvertices x 1) or (md.mesh.numberofelements x 1)") - if not min(md.qmu.partition.flatten())==0: - md.checkmessage("partition vector not indexed from 0 on") - if max(md.qmu.partition.flatten())>=md.qmu.numberofpartitions: - md.checkmessage("for qmu analysis, partitioning vector cannot go over npart, number of partition areas") + if np.size(md.qmu.vpartition) > 0: + if np.size(md.qmu.vpartition)!=md.mesh.numberofvertices: + md.checkmessage("user supplied vertex partition for qmu analysis should have size (md.mesh.numberofvertices x 1)") + if not min(md.qmu.vpartition.flatten())==0: + md.checkmessage("vertex partition vector not indexed from 0 on") + if max(md.qmu.vpartition.flatten())>=md.qmu.numberofpartitions: + md.checkmessage("for qmu analysis, vertex partitioning vector cannot go over npart, number of partition areas") + + if np.size(md.qmu.epartition) > 0: + if np.size(md.qmu.epartition) != md.mesh.numberofelements: + md.checkmessage("user supplied element partition for qmu analysis should have size (md.mesh.numberofelements x 1)") + if not min(md.qmu.epartition.flatten())==0: + md.checkmessage("elememtn partition vector not indexed from 0 on") + if max(md.qmu.epartition.flatten())>=md.qmu.numberofpartitions: + md.checkmessage("for qmu analysis, element partitioning vector cannot go over npart, number of partition areas") + + if np.size(md.qmu.vpartition) == 0 or np.any(np.isnan(md.qmu.vpartition)) or np.size(md.qmu.epartition) == 0 or np.any(np.isnan(md.qmu.epartition)): + md.checkmessage("for qmu analysis, both an element and partitioning vectors need to be supplied with no nan values! One can be defaulted to all zeros.") return md @@ -167,5 +180,6 @@ WriteData(fid,prefix,'data',False,'name','md.qmu.mass_flux_segments_present','format','Boolean'); return - WriteData(fid,prefix,'object',self,'fieldname','partition','format','DoubleMat','mattype',2) + WriteData(fid,prefix,'object',self,'fieldname','vpartition','format','DoubleMat','mattype',2) + WriteData(fid,prefix,'object',self,'fieldname','epartition','format','DoubleMat','mattype',2) WriteData(fid,prefix,'object',self,'fieldname','numberofpartitions','format','Integer') WriteData(fid,prefix,'object',self,'fieldname','numberofresponses','format','Integer') Index: /issm/trunk-jpl/src/m/partition/partitioner.m =================================================================== --- /issm/trunk-jpl/src/m/partition/partitioner.m (revision 24172) +++ /issm/trunk-jpl/src/m/partition/partitioner.m (revision 24173) @@ -9,4 +9,5 @@ % section: 1 by defaults(1=bisection, 2=quadrisection, 3=octasection) % recomputeadjacency: 'on' by default (set to 'off' to compute existing one) +% type: 'node' or 'element' partition vector (default to 'node') % Output: md.qmu.partition recover the partition vector % @@ -24,4 +25,5 @@ options=addfielddefault(options,'section',1); options=addfielddefault(options,'recomputeadjacency','on'); +options=addfielddefault(options,'type','node'); %get package: @@ -29,4 +31,5 @@ npart=getfieldvalue(options,'npart'); recomputeadjacency=getfieldvalue(options,'recomputeadjacency'); +vectortype=getfieldvalue(options,'type'); if(dimension(md.mesh)==3), @@ -53,42 +56,52 @@ if strcmpi(package,'chaco'), - % default method (from chaco.m) - method=[1 1 0 0 1 1 50 0 .001 7654321]'; - method(1)=3; % global method (3=inertial (geometric)) - method(3)=0; % vertex weights (0=off, 1=on) + if strcmpi(vectortype,'element') + error(['partitioner error message: package ' package ' does not allow element partitions.']); + else - %specify bisection - method(6)=getfieldvalue(options,'section');% ndims (1=bisection, 2=quadrisection, 3=octasection) + % default method (from chaco.m) + method=[1 1 0 0 1 1 50 0 .001 7654321]'; + method(1)=3; % global method (3=inertial (geometric)) + method(3)=0; % vertex weights (0=off, 1=on) - %are we using weights? - if strcmpi(getfieldvalue(options,'weighting'),'on'), - weights=floor(md.qmu.vertex_weight/min(md.qmu.vertex_weight)); - method(3)=1; - else - weights=[]; - end + %specify bisection + method(6)=getfieldvalue(options,'section');% ndims (1=bisection, 2=quadrisection, 3=octasection) - % partition into nparts - if isa(md.mesh,'mesh2d'), - part=Chaco(md.qmu.adjacency,weights,[],md.mesh.x, md.mesh.y,zeros(md.mesh.numberofvertices,1),method,npart,[])'+1; %index partitions from 1 up. like metis. - else - part=Chaco(md.qmu.adjacency,weights,[],md.mesh.x, md.mesh.y,md.mesh.z,method,npart,[])'+1; %index partitions from 1 up. like metis. + %are we using weights? + if strcmpi(getfieldvalue(options,'weighting'),'on'), + weights=floor(md.qmu.vertex_weight/min(md.qmu.vertex_weight)); + method(3)=1; + else + weights=[]; + end + + % partition into nparts + if isa(md.mesh,'mesh2d'), + part=Chaco(md.qmu.adjacency,weights,[],md.mesh.x, md.mesh.y,zeros(md.mesh.numberofvertices,1),method,npart,[])'+1; %index partitions from 1 up. like metis. + else + part=Chaco(md.qmu.adjacency,weights,[],md.mesh.x, md.mesh.y,md.mesh.z,method,npart,[])'+1; %index partitions from 1 up. like metis. + end + end elseif strcmpi(package,'scotch'), - %are we using weights? - if strcmpi(getfieldvalue(options,'weighting'),'on'), - weights=floor(md.qmu.vertex_weight/min(md.qmu.vertex_weight)); + if strcmpi(vectortype,'element') + error(['partitioner error message: package ' package ' does not allow element partitions.']); else - weights=[]; + %are we using weights? + if strcmpi(getfieldvalue(options,'weighting'),'on'), + weights=floor(md.qmu.vertex_weight/min(md.qmu.vertex_weight)); + else + weights=[]; + end + maptab=Scotch(md.qmu.adjacency,[],weights,[],'cmplt',[npart]); + + part=maptab(:,2)+1;%index partitions from 1 up. like metis. end - maptab=Scotch(md.qmu.adjacency,[],weights,[],'cmplt',[npart]); - - part=maptab(:,2)+1;%index partitions from 1 up. like metis. elseif strcmpi(package,'linear'), - if npart==md.mesh.numberofelements | md.qmu.numberofpartitions==md.mesh.numberofelements + if strcmpi(vectortype,'element') part=1:1:md.mesh.numberofelements; disp('Linear partitioner requesting partitions on elements'); @@ -99,5 +112,9 @@ elseif strcmpi(package,'metis'), - [element_partitioning,part]=MeshPartition(md,md.qmu.numberofpartitions); + if strcmpi(vectortype,'element') + error(['partitioner error message: package ' package ' does not allow element partitions.']); + else + [element_partitioning,part]=MeshPartition(md,md.qmu.numberofpartitions); + end else @@ -112,5 +129,9 @@ md3d.qmu.adjacency=md.qmu.adjacency; md=md3d; - part=project3d(md,'vector',part','type','node'); + if strcmpi(vectortype,'element') + part=project3d(md,'vector',part','type','element'); + else + part=project3d(md,'vector',part','type','node'); + end end @@ -118,3 +139,14 @@ part=part'; end -md.qmu.partition=part; + +if strcmpi(vectortype,'element') + md.qmu.epartition=part; + if isempty(md.qmu.vpartition) | (prod(size(md.qmu.vpartition))==1 & isnan(md.qmu.vpartition)) + md.qmu.vpartition=zeros(md.mesh.numberofvertices,1); + end +else + md.qmu.vpartition=part; + if isempty(md.qmu.epartition) | (prod(size(md.qmu.epartition))==1 & isnan(md.qmu.epartition)) + md.qmu.epartition=zeros(md.mesh.numberofelements,1); + end +end Index: /issm/trunk-jpl/src/m/partition/partitioner.py =================================================================== --- /issm/trunk-jpl/src/m/partition/partitioner.py (revision 24172) +++ /issm/trunk-jpl/src/m/partition/partitioner.py (revision 24173) @@ -21,4 +21,5 @@ section: 1 by defaults(1=bisection, 2=quadrisection, 3=octasection) recomputeadjacency: 'on' by default (set to 'off' to compute existing one) + type: 'node' or 'element' partition vector (default to 'node') Output: md.qmu.partition recover the partition vector @@ -36,4 +37,5 @@ options.addfielddefault('section',1) options.addfielddefault('recomputeadjacency','on') + options.addfielddefault('type','node') #get package: @@ -41,4 +43,8 @@ npart=options.getfieldvalue('npart') recomputeadjacency=options.getfieldvalue('recomputeadjacency') + vectortype=options.getfieldvalue('type') + + if m.strcmpi(vectortype,'element') and not m.strcmpi(package,'linear'): + raise RuntimeError('partitioner error message: package '+str(package)+' does not allow element partitions.') if(md.mesh.dimension()==3): @@ -121,8 +127,20 @@ md3d.qmu.adjacency=md.qmu.adjacency md=md3d - part=project3d(md,'vector',np.squeeze(part),'type','node') + if m.strcmpi(vectortype,'element'): + part=project3d(md,'vector',np.squeeze(part),'type','element') + else: + part=project3d(md,'vector',np.squeeze(part),'type','node') part=part.reshape(-1,1) - md.qmu.partition=part + + if m.strcmpi(vectortype,'element'): + md.qmu.epartition=part + if np.size(md.qmu.vpartition) == 0 or (np.size(md.qmu.vpartition) == 1 and np.isnan(md.qmu.vpartition)): + md.qmu.vpartition=np.zeros((md.mesh.numberofvertices,1)) + + else: + md.qmu.vpartition=part + if np.size(md.qmu.epartition) == 0 or (np.size(md.qmu.epartition) == 1 and np.isnan(md.qmu.epartition)): + md.qmu.epartition=np.zeros((md.mesh.numberofelements,1)) return md Index: /issm/trunk-jpl/src/m/qmu/importancefactors.m =================================================================== --- /issm/trunk-jpl/src/m/qmu/importancefactors.m (revision 24172) +++ /issm/trunk-jpl/src/m/qmu/importancefactors.m (revision 24173) @@ -44,7 +44,10 @@ factors=importancefactors; return; +elseif count==max(md.qmu.epartition+1) + %distribute importance factor + factors=importancefactors(md.qmu.epartition'+1); %md.qmu.partition was created to index "c" style else %distribute importance factor - factors=importancefactors(md.qmu.partition'+1); %md.qmu.partition was created to index "c" style + factors=importancefactors(md.qmu.vpartition'+1); %md.qmu.partition was created to index "c" style end Index: /issm/trunk-jpl/src/m/qmu/importancefactors.py =================================================================== --- /issm/trunk-jpl/src/m/qmu/importancefactors.py (revision 24172) +++ /issm/trunk-jpl/src/m/qmu/importancefactors.py (revision 24173) @@ -43,7 +43,11 @@ factors=importancefactors return factors + elif count==np.max(md.qmu.epartition+1): + #distribute importance factor + factors=importancefactors[(md.qmu.epartition.conj().T).flatten().astype(int)] + #md.qmu.partition was created to index "c" style else: #distribute importance factor - factors=importancefactors[(md.qmu.partition.conj().T).flatten().astype(int)] + factors=importancefactors[(md.qmu.vpartition.conj().T).flatten().astype(int)] #md.qmu.partition was created to index "c" style