C++ TPZStack::Push方法代码示例

void RefinamentoSingular(TPZAutoPointer<TPZGeoMesh> gmesh,int nref)
{
    int64_t nnodes = gmesh->NNodes();
    int64_t in;
    for (in=0; in<nnodes; in++) {
        TPZGeoNode *gno = &gmesh->NodeVec()[in];
        if (abs(gno->Coord(0))< 1.e-6 && abs(gno->Coord(1)) < 1.e-6) {
            break;
        }
    }
    if (in == nnodes) {
        DebugStop();
    }
    TPZGeoElSide gelside;
    int64_t nelem = gmesh->NElements();
    for (int64_t el = 0; el<nelem; el++) {
        TPZGeoEl *gel = gmesh->ElementVec()[el];
        int ncorner = gel->NCornerNodes();
        for (int ic=0; ic<ncorner; ic++) {
            int64_t nodeindex = gel->NodeIndex(ic);
            if (nodeindex == in) {
                gelside = TPZGeoElSide(gel, ic);
                break;
            }
        }
        if (gelside.Element()) {
            break;
        }
    }
    if (!gelside.Element()) {
        DebugStop();
    }
    for (int iref = 0; iref <nref; iref++) {
        TPZStack<TPZGeoElSide> gelstack;
        gelstack.Push(gelside);
        TPZGeoElSide neighbour = gelside.Neighbour();
        while (neighbour != gelside) {
            gelstack.Push(neighbour);
            neighbour = neighbour.Neighbour();
        }
        int64_t nstack = gelstack.size();
        for (int64_t ist=0; ist < nstack; ist++) {
            if (!gelstack[ist].Element()->HasSubElement()) {
                TPZVec<TPZGeoEl *> subel;
                gelstack[ist].Element()->Divide(subel);
            }
        }
    }
}

void TPZAnalysisError::ExpandConnected(TPZStack<TPZCompElSide> &singel){
	int64_t nelem = singel.NElements();
	TPZStack<TPZGeoElSide> gelstack;
	TPZStack<TPZCompElSide> celstack;
	int64_t iel;
	for(iel=0; iel<nelem; iel++) {
		TPZCompElSide celside = singel[iel];
		TPZGeoElSide gelside;
		gelside = celside.Reference();
		if(!gelside.Exists()) continue;
		gelstack.Resize(0);
		cout << "This part needs to be fixed\n";
		//		gelside.Element()->LowerDimensionSides(gelside.Side(),gelstack);
		while(gelstack.NElements()) {
			TPZGeoElSide gelsideloc;
			gelsideloc = gelstack.Pop();
			gelsideloc.EqualLevelCompElementList(celstack,1,0);
			while(celstack.NElements()) {
				TPZCompElSide celsideloc = celstack.Pop();
				if(! celsideloc.Exists()) continue;
				int64_t nelsing = singel.NElements();
				int64_t smel;
				for(smel=0; smel<nelsing; smel++) if(singel[smel].Element() == celsideloc.Element()) break;
				if(smel != nelsing) singel.Push(celsideloc);
			}
		}
	}
}

void TPZAdaptMesh::CreateClones(){
	// asserting references of the original meshes
    fReferenceCompMesh->Reference()->ResetReference();
    fReferenceCompMesh->LoadReferences();
    TPZGeoMesh *geomesh = fReferenceCompMesh->Reference();
    
    TPZStack<TPZGeoEl*> patch;
    
    int clid,elid;
    for (clid=0; clid<fPatchIndex.NElements()-1;clid++) {
		// making clone of the original geometric mesh, only to construct computational clone
        TPZGeoCloneMesh *geoclone = new TPZGeoCloneMesh(geomesh);
        TPZStack<TPZGeoEl*> patch;
        for (elid=fPatchIndex[clid];elid<fPatchIndex[clid+1];elid++){
            patch.Push(fPatch[elid]);
        }
        geoclone->SetElements(patch,fGeoRef[clid]);
        TPZVec<TPZGeoEl *> sub;
        //    int ngcel = geoclone->ElementVec().NElements();
        int printing = 0;
        if(printing) {
            ofstream out("testAdaptMesh.txt",ios::app);
            geoclone->Print(out);
        }
        
        TPZCompCloneMesh *clonecompmesh = new TPZCompCloneMesh(geoclone,fReferenceCompMesh);
        clonecompmesh->AutoBuild(/*fMaxP*/);
		// Computational mesh clone is stored
        fCloneMeshes.Push(clonecompmesh);    
    }
}

	void TPZPrism::LowerDimensionSides(int side,TPZStack<int> &smallsides, int DimTarget)
	{
		smallsides.Resize(0);
		int nsidecon = NContainedSides(side);
		for(int is = 0; is < nsidecon - 1; is++) {
			if (SideDimension(ContainedSideLocId(side,is)) == DimTarget) smallsides.Push(ContainedSideLocId(side,is));
		}
	}

	void TPZPrism::HigherDimensionSides(int side, TPZStack<int> &high)
	{
		if(side <0 || side >= NSides) {
			PZError << "TPZPrism::HigherDimensionSides side "<< side << endl;
		}
		int is;
		for(is=0; is<nhighdimsides[side]; is++) high.Push(highsides[side][is]);
		
	}

	void TPZRefQuad::GetSubElements(TPZGeoEl *father,int side, TPZStack<TPZGeoElSide> &subel){
		
		subel.Resize(0);
		if(side<0 || side>TPZShapeQuad::NSides || !father->HasSubElement()){
			PZError << "TPZRefQuad::GetSubelements called with error arguments\n";
			return;
		}
		int nsub = NSideSubElements(side);
		for(int i=0;i<nsub;i++)
			subel.Push(TPZGeoElSide(father->SubElement(subeldata[side][i][0]),subeldata[side][i][1]));
	}

void TPZFrontNonSym::Compress(){
	//	PrintGlobal("Before COmpress");
	//	Print("Before Compress", cout);
	TPZStack <int> from;
	int nfound;
	int i, j;
	for(i = 0; i < fFront; i++){
		if(fGlobal[i] != -1) from.Push(i);
	}
	
	/**
	 *First fLocal initialization
	 *Any needed updates is done on next loop
	 */
	nfound = from.NElements();
	for(i=0;i<nfound;i++) {
		fGlobal[i]=fGlobal[from[i]];
		//fGlobal[from[i]] = -1;
		fLocal[fGlobal[i]] = i;
	}
	for(;i<fGlobal.NElements();i++) fGlobal[i] = -1;
	
	if(nfound+fFree.NElements()!=fFront) cout << "TPZFront.Compress inconsistent data structure\n";
	
	fFront = nfound;
	fFree.Resize(0);	
	fGlobal.Resize(fFront);
	if(fData.NElements()==0) return;
	
	for(j = 0; j < nfound; j++){
		for(i = 0; i < nfound; i++){
			Element(i,j) = Element(from[i], from[j]);
			if(from[i]!=i || from[j]!=j) Element(from[i],from[j])=0.;
		}
		//		fGlobal[i] = fGlobal[from[i]];
		//		fLocal[fGlobal[i]] = i;
	}
	
	//	Print("After Compress", cout);
	//	PrintGlobal("After Compress",output);
}

//.........这里部分代码省略.........
    
    /**********************************************************************
     * Cria��o de uma malha computacional clone
     * ********************************************************************/
 	comp->GetRefPatches(patch);
	
	geomesh.ResetReference();
	TPZStack <int64_t> patchel;
	TPZStack <TPZGeoEl *> toclonegel;
	TPZStack <int64_t> patchindex;
	TPZVec<int64_t> n2elgraph;
	TPZVec<int64_t> n2elgraphid;
	TPZStack<int64_t> elgraph;
	TPZVec<int64_t> elgraphindex;
	int64_t k;
	TPZCompMesh *clonecmesh = new TPZCompMesh(&geomesh);
	cout << "Check 1: number of reference elements for patch before createcompel: " << patch.size() << endl;
    std::set<TPZGeoEl *>::iterator it;
    for (it=patch.begin(); it!=patch.end(); it++)
    {
		//patch[i]->Print(cout);
        int64_t index;
        TPZGeoEl *gel = *it;
        clonecmesh->CreateCompEl(gel, index);
        //		patch[i]->CreateCompEl(*clonecmesh,i);
	}
    //	cout << "Check 2: number of reference elements for patch after createcompel: " << patch.NElements() << endl;
	clonecmesh->CleanUpUnconnectedNodes();
    //	clonecmesh->Print(cout);
	clonecmesh->GetNodeToElGraph(n2elgraph,n2elgraphid,elgraph,elgraphindex);
	int64_t clnel = clonecmesh->NElements();
    //	cout << "Number of elements in clonemessh: " << clnel << endl;
	//o primeiro patch come�a em zero
	patchindex.Push(0);
	for (i=0; i<clnel; i++){
		//cout << endl << endl << "Evaluating patch for element: " << i << endl;
		clonecmesh->GetElementPatch(n2elgraph,n2elgraphid,elgraph,elgraphindex,i,patchel);
		cout << "Patch elements: " << patchel.NElements() << endl;
		/*for (k=0;k<patchel.NElements();k++){
         clonecmesh->ElementVec()[patchel[k]]->Reference()->Print();
         cout << endl;
         }*/
		for (j=0; j<patchel.NElements(); j++){
			//obten��o do elemento geom�trico do patch
			//cout << "Creating geometric clone elements for computational element :" << j << endl;
			TPZGeoEl *gel = clonecmesh->ElementVec()[patchel[j]]->Reference();
			//gel->Print(cout);
			//inserir todos os pais do elemento geom�trico do patch
			int64_t count = 0;
			//cout << "Inserting father element:" << "\t"; 
			while(gel){	
				TPZGeoEl *father = gel->Father();
				if (father){
					//father->Print(cout);
					gel = father;
					continue;
				}
				else toclonegel.Push(gel);
				gel = father;
				//cout <<  count << "\t";
				count ++;
			}
			//cout << endl;
		}
		int64_t sum = toclonegel.NElements()-1;
		//cout << endl << sum << endl;

//SingularElements(..)
void TPZAnalysisError::ZoomInSingularity(REAL csi, TPZCompElSide elside, REAL singularity_strength) {

	REAL hn = 1./pow(csi,1./singularity_strength);
	REAL Q=2.;
	REAL NcReal = log( 1.+(1./hn - 1.)*(Q - 1.) )/log(Q);
	int Nc = 0;
	while(REAL(Nc) < (NcReal+0.5)) Nc++;
	int minporder = 2;
	
	TPZStack<TPZCompElSide> ElToRefine;
	TPZStack<int> POrder;
	TPZStack<TPZGeoElSide> subelements;
	TPZStack<int64_t> csubindex;
	ElToRefine.Push(elside);
	POrder.Push(Nc);
	while(ElToRefine.NElements()) {
		/** Take the next element and its interpolation order from the stack*/
		TPZCompElSide curelside = ElToRefine.Pop();
		int curporder = POrder.Pop();
		if(!curelside.Exists()) continue;
		int64_t cindex = curelside.Element()->Index();
		if(cindex < 0) continue;
		
		/** Cast the element to an interpolated element if possible*/
		TPZCompEl *cel = curelside.Element();
		TPZInterpolatedElement *cintel = 0;
		cintel = dynamic_cast<TPZInterpolatedElement *> (cel);
		/** If the element is not interpolated, nothing to do */
		if(!cintel) continue;
		/** Set the interpolation order of the current element to curporder*/
		if(curporder == minporder) {
			cintel->PRefine(Nc);
			fSingular.Push(curelside);
		} else {
			cintel->PRefine(curporder);
			cintel->Divide(cindex,csubindex,1);
			/** Identify the subelements along the side and push them on the stack*/
		}
		TPZGeoElSide gelside = curelside.Reference();
		if(!gelside.Exists()) continue;
		gelside.GetSubElements2(subelements);
		int64_t ns = subelements.NElements();
		curporder--;
		int64_t is;
		for(is=0; is<ns; is++) {
			TPZGeoElSide sub = subelements[is];
			TPZCompElSide csub = sub.Reference();
			if(csub.Exists()) {
				ElToRefine.Push(csub);
				POrder.Push(curporder);
			}
		}
	}
	ExpandConnected(fSingular);
	
	/*
	 REAL H1_error,L2_error,estimate;
	 TPZBlock *flux=0;
	 int64_t nel = fElIndexes.NElements();
	 for(int64_t elloc=0;elloc<nel;elloc++) {
	 int64_t el = fElIndexes[elloc];
	 estimate = fElErrors[elloc];
	 REAL csi = estimate / fAdmissibleError;
	 REAL h = h_Parameter(intellist[el]);
	 REAL hn = h/pow(csi,1./.9);
	 REAL Nc = log( 1.+(h/hn - 1.)*(Q - 1.) )/log(Q);
	 if(hn > 1.3*h) hn = 2.0*h*hn / (h + hn);
	 REAL hsub = h;//100.0;//pode ser = h ; Cedric
	 TPZCompEl *locel = intellist[el];
	 //obter um subelemento que contem o ponto singular e tem tamanho <= hn
	 TPZAdmChunkVector<TPZCompEl *> sublist;
	 while(hsub > hn) {
	 TPZVec<int64_t> indexsubs;
	 int64_t index = locel->Index();
	 locel->Divide(index,indexsubs,1);
	 int64_t nsub = indexsubs.NElements();
	 TPZAdmChunkVector<TPZCompEl *> listsub(0);
	 for(int64_t k=0;k<nsub;k++) {
	 index = listsub.AllocateNewElement();
	 listsub[index] = Mesh()->ElementVec()[indexsubs[k]];
	 }
	 //existe um unico filho que contem o ponto singular
	 SingularElement(point,listsub,sublist);
	 hsub = h_Parameter(sublist[0]);
	 }
	 TPZInterpolatedElement *intel = (TPZInterpolatedElement *) locel;
	 intel->PRefine(Nc+1);
	 indexlist.Push(intel->Index());
	 //os elemento viz devem ter ordens menores a cel quanto mais longe de point
	 TPZInterpolatedElement *neighkeep,*neigh;
	 //feito s�para o caso 1d , extender para o caso geral
	 int dim = intel->Dimension();
	 if(dim != 1) {
	 cout << "TPZAnalysisError::Step3 not dimension implemented , dimension = " << intellist[el]->Dimension() << endl;
	 return ;//exit(1);
	 }
	 for(int side=0;side<2;side++) {
	 int ly = 1;
	 TPZGeoElSide neighside = intel->Reference()->Neighbour(side);
//.........这里部分代码省略.........