From: ageay Date: Wed, 28 Nov 2012 16:09:45 +0000 (+0000) Subject: warning hunting + tests on gauss pt RW on profiles with different discretizations... X-Git-Tag: V6_6_0~4 X-Git-Url: http://git.salome-platform.org/gitweb/?a=commitdiff_plain;h=1473ae61df8648d7cd850088b229b97ffa9e3248;p=tools%2Fmedcoupling.git warning hunting + tests on gauss pt RW on profiles with different discretizations per geometric type. --- diff --git a/src/MEDLoader/MEDFileData.cxx b/src/MEDLoader/MEDFileData.cxx index b18c3cb42..61f1429f3 100644 --- a/src/MEDLoader/MEDFileData.cxx +++ b/src/MEDLoader/MEDFileData.cxx @@ -129,7 +129,6 @@ bool MEDFileData::unPolyzeMeshes() throw(INTERP_KERNEL::Exception) MEDFileMeshes *ms=_meshes; if(!ms) return false; - bool ret=false; std::vector< MEDFileMesh * > meshesImpacted; std::vector< DataArrayInt * > renumParamsOfMeshImpacted;//same size as meshesImpacted std::vector< std::vector > oldCodeOfMeshImpacted,newCodeOfMeshImpacted;//same size as meshesImpacted diff --git a/src/MEDLoader/MEDFileField.cxx b/src/MEDLoader/MEDFileField.cxx index e24fb86d2..70d126a7d 100644 --- a/src/MEDLoader/MEDFileField.cxx +++ b/src/MEDLoader/MEDFileField.cxx @@ -254,18 +254,32 @@ void MEDFileFieldPerMeshPerTypePerDisc::assignFieldNoProfile(int& start, int off } /*! - * Leaf method of field with profile assignement. - * @param pflName input containing name of profile if any. 0 if no profile. - * @param multiTypePfl input containing the profile array \b including \b all \b types. This array is usefull only for GAUSS_NE. - * @param idsInPfl input containing the ids in the profile 'multiTypePfl' concerning the current geo type. + * Leaf method of field with profile assignement. This method is the most general one. No optimization is done here. + * \param [in] pflName input containing name of profile if any. 0 if no profile (except for GAUSS_PT where a no profile can hide a profile when splitted by loc_id). + * \param [in] multiTypePfl is the end user profile specified in high level API + * \param [in] idsInPfl is the selection into the \a multiTypePfl whole profile that corresponds to the current geometric type. + * \param [in] locIds is the profile needed to be created for MED file format. It can be null if all cells of current geometric type are fetched in \a multiTypePfl. + * \b WARNING if not null the MED file profile can be subdivided again in case of Gauss points. + * \param [in] mesh is the mesh coming from the MEDFileMesh instance in correspondance with the MEDFileField. The mesh inside the \a field is simply ignored. */ -void MEDFileFieldPerMeshPerTypePerDisc::assignFieldProfile(int& start, const char *pflName, const DataArrayInt *multiTypePfl, const DataArrayInt *idsInPfl, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) +void MEDFileFieldPerMeshPerTypePerDisc::assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const DataArrayInt *idsInPfl, DataArrayInt *locIds, int nbOfEltsInWholeMesh, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) { - if(pflName) - _profile=pflName; - else - _profile.clear(); + _profile.clear(); _type=field->getTypeOfField(); + std::string pflName(multiTypePfl->getName()); + std::ostringstream oss; oss << pflName; + if(_type!=ON_NODES) { const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(getGeoType()); oss << "_" << cm.getRepr(); } else { oss << "_NODE"; } + if(locIds) + { + if(pflName.empty()) + throw INTERP_KERNEL::Exception("MEDFileFieldPerMeshPerTypePerDisc::assignFieldProfile : existing profile with empty name !"); + if(_type!=ON_GAUSS_PT) + { + locIds->setName(oss.str().c_str()); + glob.appendProfile(locIds); + _profile=oss.str(); + } + } const DataArrayDouble *da=field->getArray(); _start=start; switch(_type) @@ -288,7 +302,7 @@ void MEDFileFieldPerMeshPerTypePerDisc::assignFieldProfile(int& start, const cha { MEDCouplingAutoRefCountObjectPtr arr=field->getDiscretization()->getOffsetArr(mesh); MEDCouplingAutoRefCountObjectPtr arr2=arr->deltaShiftIndex(); - MEDCouplingAutoRefCountObjectPtr arr3=arr2->selectByTupleId(multiTypePfl->getConstPointer(),multiTypePfl->getConstPointer()+multiTypePfl->getNumberOfTuples()); + MEDCouplingAutoRefCountObjectPtr arr3=arr2->selectByTupleId(multiTypePfl->begin(),multiTypePfl->end()); arr3->computeOffsets2(); MEDCouplingAutoRefCountObjectPtr tmp=idsInPfl->buildExplicitArrByRanges(arr3); int trueNval=tmp->getNumberOfTuples(); @@ -299,7 +313,52 @@ void MEDFileFieldPerMeshPerTypePerDisc::assignFieldProfile(int& start, const cha } case ON_GAUSS_PT: { - throw INTERP_KERNEL::Exception("MEDFileFieldPerMeshPerTypePerDisc::assignFieldProfile : not implemented yet for profiles on gauss points !"); + const MEDCouplingFieldDiscretizationGauss *disc2=dynamic_cast(field->getDiscretization()); + if(!disc2) + throw INTERP_KERNEL::Exception("addNewEntryIfNecessaryGauss : invalid call to this method ! Internal Error !"); + const DataArrayInt *da1=disc2->getArrayOfDiscIds(); + const MEDCouplingGaussLocalization& gsLoc=field->getGaussLocalization(_loc_id); + MEDCouplingAutoRefCountObjectPtr da2=da1->selectByTupleId(idsInPfl->begin(),idsInPfl->end()); + MEDCouplingAutoRefCountObjectPtr da3=da2->getIdsEqual(_loc_id); + MEDCouplingAutoRefCountObjectPtr da4=idsInPfl->selectByTupleId(da3->begin(),da3->end()); + // + MEDCouplingAutoRefCountObjectPtr mesh2=mesh->buildPart(multiTypePfl->begin(),multiTypePfl->end()); + MEDCouplingAutoRefCountObjectPtr arr=disc2->getOffsetArr(mesh2); + // + MEDCouplingAutoRefCountObjectPtr tmp=DataArrayInt::New(); + int trueNval=0; + for(const int *pt=da4->begin();pt!=da4->end();pt++) + trueNval+=arr->getIJ(*pt+1,0)-arr->getIJ(*pt,0); + tmp->alloc(trueNval,1); + int *tmpPtr=tmp->getPointer(); + for(const int *pt=da4->begin();pt!=da4->end();pt++) + for(int j=arr->getIJ(*pt,0);jgetIJ(*pt+1,0);j++) + *tmpPtr++=j; + // + _nval=da4->getNumberOfTuples(); + getArray()->setContigPartOfSelectedValues(_start,da,tmp); + _end=_start+trueNval; + oss << "_loc_" << _loc_id; + if(locIds) + { + MEDCouplingAutoRefCountObjectPtr da5=locIds->selectByTupleId(da3->begin(),da3->end()); + da5->setName(oss.str().c_str()); + glob.appendProfile(da5); + _profile=oss.str(); + } + else + { + if(da3->getNumberOfTuples()!=nbOfEltsInWholeMesh || !da3->isIdentity()) + { + da3->setName(oss.str().c_str()); + glob.appendProfile(da3); + _profile=oss.str(); + } + } + std::ostringstream oss2; oss2 << "Loc_" << getName() << "_" << INTERP_KERNEL::CellModel::GetCellModel(getGeoType()).getRepr() << "_" << _loc_id; + _localization=oss2.str(); + glob.appendLoc(_localization.c_str(),getGeoType(),gsLoc.getRefCoords(),gsLoc.getGaussCoords(),gsLoc.getWeights()); + break; } default: throw INTERP_KERNEL::Exception("MEDFileFieldPerMeshPerTypePerDisc::assignFieldProfile : not implemented yet for such discretization type of field !"); @@ -620,10 +679,10 @@ int MEDFileFieldPerMeshPerTypePerDisc::fillEltIdsFromCode(int offset, const std: _loc_id=offset; std::ostringstream oss; std::size_t nbOfType=codeOfMesh.size()/3; - std::size_t found=-1; + int found=-1; for(std::size_t i=0;igetNumberOfTuples() << " whereas the number of ids is set to " << _nval << " for this geometric type !"; throw INTERP_KERNEL::Exception(oss.str().c_str()); } - int offset=codeOfMesh[3*found+2]; + int offset2=codeOfMesh[3*found+2]; for(const int *pflId=pfl->begin();pflId!=pfl->end();pflId++) { if(*pflId_loc_id)==geoType && (*it)->_nval==nbMeshEntities) { if(!isPfl) - if((*it)->_profile.empty()) - break; - else - if(!(*it)->_profile.empty()) - { - const DataArrayInt *pfl=glob.getProfile((*it)->_profile.c_str()); - if(pfl->isEqualWithoutConsideringStr(*idsOfMeshElt)) - break; - } + { + if((*it)->_profile.empty()) + break; + else + if(!(*it)->_profile.empty()) + { + const DataArrayInt *pfl=glob.getProfile((*it)->_profile.c_str()); + if(pfl->isEqualWithoutConsideringStr(*idsOfMeshElt)) + break; + } + } } } if(it==entriesOnSameDisc.end()) @@ -865,28 +926,20 @@ void MEDFileFieldPerMeshPerType::assignFieldNoProfile(int& start, int offset, in _field_pm_pt_pd[*it]->assignFieldNoProfile(start,offset,nbOfCells,field,glob); } -void MEDFileFieldPerMeshPerType::assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const DataArrayInt *idsInPfl, DataArrayInt *locIds, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) +/*! + * This method is the most general one. No optimization is done here. + * \param [in] multiTypePfl is the end user profile specified in high level API + * \param [in] idsInPfl is the selection into the \a multiTypePfl whole profile that corresponds to the current geometric type. + * \param [in] locIds is the profile needed to be created for MED file format. It can be null if all cells of current geometric type are fetched in \a multiTypePfl. + * \b WARNING if not null the MED file profile can be subdivided again in case of Gauss points. + * \param [in] nbOfEltsInWholeMesh nb of elts of type \a this->_geo_type in \b WHOLE mesh + * \param [in] mesh is the mesh coming from the MEDFileMesh instance in correspondance with the MEDFileField. The mesh inside the \a field is simply ignored. + */ +void MEDFileFieldPerMeshPerType::assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const DataArrayInt *idsInPfl, DataArrayInt *locIds, int nbOfEltsInWholeMesh, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) { std::vector pos=addNewEntryIfNecessary(field,idsInPfl); - if(locIds) - { - // - std::string pflName(locIds->getName()); - if(pflName.empty()) - throw INTERP_KERNEL::Exception("MEDFileFieldPerMeshPerType::assignFieldProfile : existing profile with empty name !"); - const INTERP_KERNEL::CellModel& cm=INTERP_KERNEL::CellModel::GetCellModel(_geo_type); - std::ostringstream oss; oss << pflName << "_" << cm.getRepr(); - locIds->setName(oss.str().c_str()); - glob.appendProfile(locIds); - // - for(std::vector::const_iterator it=pos.begin();it!=pos.end();it++) - _field_pm_pt_pd[*it]->assignFieldProfile(start,oss.str().c_str(),multiTypePfl,idsInPfl,field,mesh,glob); - } - else - { - for(std::vector::const_iterator it=pos.begin();it!=pos.end();it++) - _field_pm_pt_pd[*it]->assignFieldProfile(start,0,multiTypePfl,idsInPfl,field,mesh,glob); - } + for(std::vector::const_iterator it=pos.begin();it!=pos.end();it++) + _field_pm_pt_pd[*it]->assignFieldProfile(start,multiTypePfl,idsInPfl,locIds,nbOfEltsInWholeMesh,field,mesh,glob); } void MEDFileFieldPerMeshPerType::assignNodeFieldNoProfile(int& start, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) @@ -898,17 +951,11 @@ void MEDFileFieldPerMeshPerType::assignNodeFieldNoProfile(int& start, const MEDC void MEDFileFieldPerMeshPerType::assignNodeFieldProfile(int& start, const DataArrayInt *pfl, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) { - std::string pflName(pfl->getName()); - if(pflName.empty()) - throw INTERP_KERNEL::Exception("MEDFileFieldPerMeshPerType::assignNodeFieldProfile : existing profile with empty name !"); - std::ostringstream oss; oss << pflName << "_NODE"; MEDCouplingAutoRefCountObjectPtr pfl2=pfl->deepCpy(); - pfl2->setName(oss.str().c_str()); - glob.appendProfile(pfl2); // _field_pm_pt_pd.resize(1); _field_pm_pt_pd[0]=MEDFileFieldPerMeshPerTypePerDisc::New(this,ON_NODES,-3); - _field_pm_pt_pd[0]->assignFieldProfile(start,oss.str().c_str(),pfl,pfl2,field,0,glob);//mesh is not requested so 0 is send. + _field_pm_pt_pd[0]->assignFieldProfile(start,pfl,pfl2,pfl2,-1,field,0,glob);//mesh is not requested so 0 is send. } std::vector MEDFileFieldPerMeshPerType::addNewEntryIfNecessary(const MEDCouplingFieldDouble *field, int offset, int nbOfCells) throw(INTERP_KERNEL::Exception) @@ -1043,7 +1090,7 @@ std::vector MEDFileFieldPerMeshPerType::addNewEntryIfNecessaryGauss(const M if(!disc2) throw INTERP_KERNEL::Exception("addNewEntryIfNecessaryGauss : invalid call to this method ! Internal Error !"); const DataArrayInt *da=disc2->getArrayOfDiscIds(); - MEDCouplingAutoRefCountObjectPtr da2=da->selectByTupleId(subCells->getConstPointer(),subCells->getConstPointer()+subCells->getNumberOfTuples()); + MEDCouplingAutoRefCountObjectPtr da2=da->selectByTupleIdSafe(subCells->getConstPointer(),subCells->getConstPointer()+subCells->getNumberOfTuples()); std::set retTmp=da2->getDifferentValues(); if(retTmp.find(-1)!=retTmp.end()) throw INTERP_KERNEL::Exception("addNewEntryIfNecessaryGauss : some cells have no dicretization description !"); @@ -1403,24 +1450,6 @@ void MEDFileFieldPerMesh::copyTinyInfoFrom(const MEDCouplingMesh *mesh) throw(IN mesh->getTime(_mesh_iteration,_mesh_order); } -void MEDFileFieldPerMesh::assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const std::vector& code, const std::vector& idsInPflPerType, const std::vector& idsPerType, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) -{ - int nbOfTypes=code.size()/3; - bool isProfile=false; - for(int i=0;i& code, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) { int nbOfTypes=code.size()/3; @@ -1437,8 +1466,14 @@ void MEDFileFieldPerMesh::assignFieldNoProfileNoRenum(int& start, const std::vec /*! * This method is the most general one. No optimization is done here. + * \param [in] multiTypePfl is the end user profile specified in high level API + * \param [in] code is the code of \a mesh[multiTypePfl] mesh. It is of size of number of different geometric types into \a mesh[multiTypePfl]. + * \param [in] code2 is the code of the \b WHOLE mesh on the same level. So all types in \a code are in \a code2. + * \param [in] idsInPflPerType is the selection into the \a multiTypePfl whole profile that corresponds to the given geometric type. This vector is always 3 times smaller than \a code. + * \param [in] idsPerType is a vector containing the profiles needed to be created for MED file format. \b WARNING these processed MED file profiles can be subdivided again in case of Gauss points. + * \param [in] mesh is the mesh coming from the MEDFileMesh instance in correspondance with the MEDFileField. The mesh inside the \a field is simply ignored. */ -void MEDFileFieldPerMesh::assignFieldProfileGeneral(int& start, const DataArrayInt *multiTypePfl, const std::vector& code, const std::vector& idsInPflPerType, const std::vector& idsPerType, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) +void MEDFileFieldPerMesh::assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const std::vector& code, const std::vector& code2, const std::vector& idsInPflPerType, const std::vector& idsPerType, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception) { int nbOfTypes=code.size()/3; for(int i=0;iassignFieldProfile(start,multiTypePfl,idsInPflPerType[i],pfl,field,mesh,glob); + int nbOfTupes2=code2.size()/3; + int found=0; + for(;foundassignFieldProfile(start,multiTypePfl,idsInPflPerType[i],pfl,code2[3*found+1],field,mesh,glob); } } @@ -1831,6 +1873,9 @@ void MEDFileFieldPerMesh::changeLocsRefsNamesGen(const std::vector< std::pairchangeLocsRefsNamesGen(mapOfModif); } +/*! + * \param [in] mesh is the whole mesh + */ MEDCouplingFieldDouble *MEDFileFieldPerMesh::getFieldOnMeshAtLevel(TypeOfField type, const MEDFileFieldGlobsReal *glob, const MEDCouplingMesh *mesh, bool& isPfl) const throw(INTERP_KERNEL::Exception) { if(_field_pm_pt.empty()) @@ -2006,7 +2051,8 @@ int MEDFileFieldPerMesh::addNewEntryIfNecessary(INTERP_KERNEL::NormalizedCellTyp } /*! - * 'dads' and 'locs' input parameters have the same number of elements. + * 'dads' and 'locs' input parameters have the same number of elements + * \param [in] mesh is \b NOT the global mesh, but the possibly reduced mesh. \a mesh parameter will be directly aggregated in the returned field */ MEDCouplingFieldDouble *MEDFileFieldPerMesh::finishField(TypeOfField type, const MEDFileFieldGlobsReal *glob, const std::vector< std::pair >& dads, const std::vector& locs, @@ -2047,6 +2093,7 @@ MEDCouplingFieldDouble *MEDFileFieldPerMesh::finishField(TypeOfField type, const * 'dads', 'locs' and 'geoTypes' input parameters have the same number of elements. * No check of this is performed. 'da' array contains an array in old2New style to be applyied to mesh to obtain the right support. * The order of cells in the returned field is those imposed by the profile. + * \param [in] mesh is the global mesh. */ MEDCouplingFieldDouble *MEDFileFieldPerMesh::finishField2(TypeOfField type, const MEDFileFieldGlobsReal *glob, const std::vector >& dads, const std::vector& locs, @@ -2059,11 +2106,10 @@ MEDCouplingFieldDouble *MEDFileFieldPerMesh::finishField2(TypeOfField type, cons if(nbOfTuples==mesh->getNumberOfCells()) return finishField(type,glob,dads,locs,mesh,isPfl); } - MEDCouplingAutoRefCountObjectPtr ret=finishField(type,glob,dads,locs,mesh,isPfl); - isPfl=true; MEDCouplingAutoRefCountObjectPtr m2=mesh->buildPart(da->getConstPointer(),da->getConstPointer()+da->getNbOfElems()); m2->setName(mesh->getName()); - ret->setMesh(m2); + MEDCouplingAutoRefCountObjectPtr ret=finishField(type,glob,dads,locs,m2,isPfl); + isPfl=true; ret->incrRef(); return ret; } @@ -3420,18 +3466,14 @@ void MEDFileField1TSWithoutSDA::setFieldNoProfileSBT(const MEDCouplingFieldDoubl TypeOfField type=field->getTypeOfField(); std::vector dummy; int start=copyTinyInfoFrom(field); + int pos=addNewEntryIfNecessary(mesh); if(type!=ON_NODES) { std::vector code=MEDFileField1TSWithoutSDA::CheckSBTMesh(mesh); - // - int pos=addNewEntryIfNecessary(mesh); - _field_per_mesh[pos]->assignFieldProfile(start,0,code,dummy,dummy,field,0,glob);//mesh is set to 0 because no external mesh is needed to be sent because no profile. + _field_per_mesh[pos]->assignFieldNoProfileNoRenum(start,code,field,glob); } else - { - int pos=addNewEntryIfNecessary(mesh); - _field_per_mesh[pos]->assignNodeFieldNoProfile(start,field,glob); - } + _field_per_mesh[pos]->assignNodeFieldNoProfile(start,field,glob); } /*! @@ -3443,11 +3485,12 @@ void MEDFileField1TSWithoutSDA::setFieldProfile(const MEDCouplingFieldDouble *fi int start=copyTinyInfoFrom(field); std::vector idsInPflPerType; std::vector idsPerType; - std::vector code; + std::vector code,code2; MEDCouplingAutoRefCountObjectPtr m=mesh->getGenMeshAtLevel(meshDimRelToMax); if(type!=ON_NODES) { m->splitProfilePerType(profile,code,idsInPflPerType,idsPerType); + code2=m->getDistributionOfTypes(); // std::vector< MEDCouplingAutoRefCountObjectPtr > idsInPflPerType2(idsInPflPerType.size()); for(std::size_t i=0;iassignFieldProfile(start,profile,code,idsInPflPerType,idsPerType,field,m,glob); + _field_per_mesh[pos]->assignFieldProfile(start,profile,code,code2,idsInPflPerType,idsPerType,field,m,glob); } else { @@ -3498,6 +3541,9 @@ MEDCouplingFieldDouble *MEDFileField1TSWithoutSDA::getFieldAtTopLevel(TypeOfFiel return MEDFileField1TSWithoutSDA::getFieldOnMeshAtLevel(type,meshDimRelToMax,renumPol,glob,mm); } +/*! + * \param [in] mesh is the whole mesh. + */ MEDCouplingFieldDouble *MEDFileField1TSWithoutSDA::getFieldOnMeshAtLevel(TypeOfField type, int renumPol, const MEDFileFieldGlobsReal *glob, const MEDCouplingMesh *mesh, const DataArrayInt *cellRenum, const DataArrayInt *nodeRenum) const throw(INTERP_KERNEL::Exception) { static const char msg1[]="MEDFileField1TSWithoutSDA::getFieldOnMeshAtLevel : request for a renumbered field following mesh numbering whereas it is a profile field !"; @@ -3597,7 +3643,7 @@ DataArrayDouble *MEDFileField1TSWithoutSDA::getUndergroundDataArrayExt(std::vect } MEDFileField1TSWithoutSDA::MEDFileField1TSWithoutSDA(const char *fieldName, int csit, int fieldtype, int iteration, int order, - const std::vector& infos):_csit(csit),_field_type(fieldtype),_iteration(iteration),_order(order) + const std::vector& infos):_iteration(iteration),_order(order),_csit(csit),_field_type(fieldtype) { DataArrayDouble *arr=getOrCreateAndGetArray(); arr->setName(fieldName); @@ -3747,7 +3793,7 @@ void MEDFileField1TS::write(const char *fileName, int mode) const throw(INTERP_K } MEDFileField1TS::MEDFileField1TS(const char *fileName, const char *fieldName, int iteration, int order) throw(INTERP_KERNEL::Exception) -try:_content(MEDFileField1TSWithoutSDA::New(fieldName,-1,-1,iteration,order,std::vector())),MEDFileFieldGlobsReal(fileName) +try:MEDFileFieldGlobsReal(fileName),_content(MEDFileField1TSWithoutSDA::New(fieldName,-1,-1,iteration,order,std::vector())) { MEDFileUtilities::CheckFileForRead(fileName); MEDFileUtilities::AutoFid fid=MEDfileOpen(fileName,MED_ACC_RDONLY); @@ -4218,7 +4264,6 @@ MEDFileFieldMultiTSWithoutSDA::MEDFileFieldMultiTSWithoutSDA(med_idt fid, int fi try:_name("") { med_field_type typcha; - int nbstep2=-1; // int ncomp=MEDfieldnComponent(fid,fieldId+1); INTERP_KERNEL::AutoPtr comp=MEDLoaderBase::buildEmptyString(ncomp*MED_SNAME_SIZE); @@ -4539,12 +4584,12 @@ int MEDFileFieldMultiTSWithoutSDA::getPosOfTimeStep(int iteration, int order) co const MEDFileField1TSWithoutSDA *tmp(*it); if(tmp) { - int it,ord; - tmp->getTime(it,ord); - if(it==iteration && order==ord) + int it2,ord; + tmp->getTime(it2,ord); + if(it2==iteration && order==ord) return ret; else - oss << "(" << it << "," << ord << "), "; + oss << "(" << it2 << "," << ord << "), "; } } throw INTERP_KERNEL::Exception(oss.str().c_str()); @@ -4560,8 +4605,8 @@ int MEDFileFieldMultiTSWithoutSDA::getPosGivenTime(double time, double eps) cons const MEDFileField1TSWithoutSDA *tmp(*it); if(tmp) { - int it,ord; - double ti=tmp->getTime(it,ord); + int it2,ord; + double ti=tmp->getTime(it2,ord); if(fabs(time-ti)& dad, std::string& pfl, std::string& loc) const throw(INTERP_KERNEL::Exception); void writeLL(med_idt fid) const throw(INTERP_KERNEL::Exception); @@ -165,7 +165,7 @@ namespace ParaMEDMEM static MEDFileFieldPerMeshPerType *New(MEDFileFieldPerMesh *fath, INTERP_KERNEL::NormalizedCellType geoType) throw(INTERP_KERNEL::Exception); static MEDFileFieldPerMeshPerType *NewOnRead(med_idt fid, MEDFileFieldPerMesh *fath, TypeOfField type, INTERP_KERNEL::NormalizedCellType geoType) throw(INTERP_KERNEL::Exception); void assignFieldNoProfile(int& start, int offset, int nbOfCells, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); - void assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const DataArrayInt *idsInPfl, DataArrayInt *locIds, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); + void assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const DataArrayInt *idsInPfl, DataArrayInt *locIds, int nbOfEltsInWholeMesh, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); void assignNodeFieldNoProfile(int& start, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); void assignNodeFieldProfile(int& start, const DataArrayInt *pfl, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); const MEDFileFieldPerMesh *getFather() const; @@ -219,8 +219,7 @@ namespace ParaMEDMEM static MEDFileFieldPerMesh *NewOnRead(med_idt fid, MEDFileField1TSWithoutSDA *fath, int meshCsit, int meshIteration, int meshOrder) throw(INTERP_KERNEL::Exception); void simpleRepr(int bkOffset,std::ostream& oss, int id) const; void copyTinyInfoFrom(const MEDCouplingMesh *mesh) throw(INTERP_KERNEL::Exception); - void assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const std::vector& code, const std::vector& idsInPflPerType, const std::vector& idsPerType, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); - void assignFieldProfileGeneral(int& start, const DataArrayInt *multiTypePfl, const std::vector& code, const std::vector& idsInPflPerType, const std::vector& idsPerType, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); + void assignFieldProfile(int& start, const DataArrayInt *multiTypePfl, const std::vector& code, const std::vector& code2, const std::vector& idsInPflPerType, const std::vector& idsPerType, const MEDCouplingFieldDouble *field, const MEDCouplingMesh *mesh, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); void assignFieldNoProfileNoRenum(int& start, const std::vector& code, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); void assignNodeFieldNoProfile(int& start, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); void assignNodeFieldProfile(int& start, const DataArrayInt *pfl, const MEDCouplingFieldDouble *field, MEDFileFieldGlobsReal& glob) throw(INTERP_KERNEL::Exception); diff --git a/src/MEDLoader/MEDFileMesh.cxx b/src/MEDLoader/MEDFileMesh.cxx index 89620ee73..5c5ccd8f8 100644 --- a/src/MEDLoader/MEDFileMesh.cxx +++ b/src/MEDLoader/MEDFileMesh.cxx @@ -217,8 +217,8 @@ std::vector MEDFileMesh::getFamiliesOnGroups(const std::vector grps=getGroupsNames(); oss << "\" !\nAvailable groups are :"; - std::copy(grps.begin(),grps.end(),std::ostream_iterator(oss," ")); + std::vector grps2=getGroupsNames(); oss << "\" !\nAvailable groups are :"; + std::copy(grps2.begin(),grps2.end(),std::ostream_iterator(oss," ")); throw INTERP_KERNEL::Exception(oss.str().c_str()); } fams.insert((*it2).second.begin(),(*it2).second.end()); diff --git a/src/MEDLoader/MEDFileMeshLL.cxx b/src/MEDLoader/MEDFileMeshLL.cxx index 2c3c77d9b..cfb86ad46 100644 --- a/src/MEDLoader/MEDFileMeshLL.cxx +++ b/src/MEDLoader/MEDFileMeshLL.cxx @@ -652,7 +652,6 @@ void MEDFileUMeshSplitL1::eraseFamilyField() void MEDFileUMeshSplitL1::setGroupsFromScratch(const std::vector& ms, std::map& familyIds, std::map >& groups) throw(INTERP_KERNEL::Exception) { - int sz=ms.size(); std::vector< DataArrayInt * > corr; _m=MEDCouplingUMesh::FuseUMeshesOnSameCoords(ms,0,corr); std::vector< MEDCouplingAutoRefCountObjectPtr > corrMSafe(corr.begin(),corr.end()); diff --git a/src/MEDLoader/MEDLoader.cxx b/src/MEDLoader/MEDLoader.cxx index 60e8773ec..842695a46 100644 --- a/src/MEDLoader/MEDLoader.cxx +++ b/src/MEDLoader/MEDLoader.cxx @@ -261,11 +261,11 @@ void MEDLoader::MEDConnOfOneElemType::releaseArray() delete [] _global; } -MEDLoader::MEDFieldDoublePerCellType::MEDFieldDoublePerCellType(INTERP_KERNEL::NormalizedCellType type, double *values, int ncomp, int ntuple, - const int *cellIdPerType, const char *locName):_ntuple(ntuple),_ncomp(ncomp),_values(values),_type(type) +MEDLoader::MEDFieldDoublePerCellType::MEDFieldDoublePerCellType(INTERP_KERNEL::NormalizedCellType type, double *values, int ncomp, int nGeoElt, int nbi, + const int *cellIdPerType, const char *locName):_ngeo_elt(nGeoElt),_nbi(nbi),_ncomp(ncomp),_values(values),_type(type) { if(cellIdPerType) - _cell_id_per_type.insert(_cell_id_per_type.end(),cellIdPerType,cellIdPerType+ntuple); + _cell_id_per_type.insert(_cell_id_per_type.end(),cellIdPerType,cellIdPerType+nGeoElt); if(locName) _loc_name=locName; } @@ -308,6 +308,7 @@ void MEDLoaderNS::fillGaussDataOnField(const char *fileName, const std::list::const_iterator iter=data.begin();iter!=data.end();iter++) { const std::string& loc=(*iter).getLocName(); @@ -329,7 +330,16 @@ void MEDLoaderNS::fillGaussDataOnField(const char *fileName, const std::list refcoo(nbPtPerCell*dim),gscoo(nbOfGaussPt*dim),w(nbOfGaussPt); MEDlocalizationRd(fid,(*iter).getLocName().c_str(),MED_FULL_INTERLACE,&refcoo[0],&gscoo[0],&w[0]); - f->setGaussLocalizationOnType((*iter).getType(),refcoo,gscoo,w); + if((*iter).getCellIdPerType().empty()) + f->setGaussLocalizationOnType((*iter).getType(),refcoo,gscoo,w); + else + { + MEDCouplingAutoRefCountObjectPtr pfl=DataArrayInt::New(); + pfl->alloc((*iter).getCellIdPerType().size(),1); + pfl->iota(offset); + f->setGaussLocalizationOnCells(pfl->begin(),pfl->end(),refcoo,gscoo,w); + } + offset+=(*iter).getNbOfGeoElt(); } MEDfileClose(fid); } @@ -1106,27 +1116,33 @@ void MEDLoaderNS::readFieldDoubleDataInMedFile(const char *fileName, const char { if(nbPdt>0) { - int profilesize,nbi; - int nval=MEDfieldnValueWithProfile(fid,fieldName,numdt,numo,tabEnt[typeOfOutField],tabType[typeOfOutField][j],1,MED_COMPACT_PFLMODE,pflname,&profilesize,locname,&nbi); - if(nval>0) + INTERP_KERNEL::AutoPtr pflDummy=MEDLoaderBase::buildEmptyString(MED_NAME_SIZE); + INTERP_KERNEL::AutoPtr locDummy=MEDLoaderBase::buildEmptyString(MED_NAME_SIZE); + int nbProfiles=MEDfieldnProfile(fid,fieldName,numdt,numo,tabEnt[typeOfOutField],tabType[typeOfOutField][j],pflDummy,locDummy); + for(int kk=0;kk0) { - delete [] valr; - continue; + double *valr=new double[ncomp*nval*nbi]; + MEDfieldValueWithProfileRd(fid,fieldName,iteration,order,tabEnt[typeOfOutField],tabType[typeOfOutField][j],MED_COMPACT_PFLMODE, + pflname,MED_FULL_INTERLACE,MED_ALL_CONSTITUENT,(unsigned char*)valr); + std::string tmp(locname); + if((locname[0]!='\0' && (typeOfOutField!=ON_GAUSS_PT)) + || (locname[0]=='\0' && typeOfOutField==ON_GAUSS_PT)) + { + delete [] valr; + continue; + } + INTERP_KERNEL::AutoPtr pfl=0; + if(pflname[0]!='\0') + { + pfl=new int[nval]; + MEDprofileRd(fid,pflname,pfl); + } + field.push_back(MEDLoader::MEDFieldDoublePerCellType(typmai2[j],valr,ncomp,nval,nbi,pfl,locname)); } - INTERP_KERNEL::AutoPtr pfl=0; - if(pflname[0]!='\0') - { - pfl=new int[nval]; - MEDprofileRd(fid,pflname,pfl); - } - field.push_back(MEDLoader::MEDFieldDoublePerCellType(typmai2[j],valr,ncomp,nval*nbi,pfl,locname)); } } } @@ -1884,32 +1900,26 @@ ParaMEDMEM::MEDCouplingFieldDouble *MEDLoaderNS::readFieldDoubleLev2(const char throw INTERP_KERNEL::Exception(oss.str().c_str()); } //for profiles - ParaMEDMEM::MEDCouplingUMesh *newMesh=0; + MEDCouplingAutoRefCountObjectPtr newMesh; std::string mName(mesh->getName()); - for(std::list::const_iterator iter=fieldPerCellType.begin();iter!=fieldPerCellType.end();iter++) + if(typeOfOutField==ON_NODES) { - const std::vector& cellIds=(*iter).getCellIdPerType(); - if(!cellIds.empty()) + for(std::list::const_iterator iter=fieldPerCellType.begin();iter!=fieldPerCellType.end();iter++) { - std::vector ci(cellIds.size()); - std::transform(cellIds.begin(),cellIds.end(),ci.begin(),std::bind2nd(std::plus(),-1)); - ParaMEDMEM::MEDCouplingUMesh *mesh2=0; - if(typeOfOutField==ON_CELLS) + const std::vector& cellIds=(*iter).getCellIdPerType(); + if(!cellIds.empty()) { - if(newMesh) - mesh2=newMesh->keepSpecifiedCells((*iter).getType(),&ci[0],&ci[0]+ci.size()); - else - mesh2=mesh->keepSpecifiedCells((*iter).getType(),&ci[0],&ci[0]+ci.size()); - } - else if(typeOfOutField==ON_NODES) - { - DataArrayInt *da=0,*da2=0; - if(newMesh) + std::vector ci(cellIds.size()); + std::transform(cellIds.begin(),cellIds.end(),ci.begin(),std::bind2nd(std::plus(),-1)); + MEDCouplingAutoRefCountObjectPtr mesh2; + MEDCouplingAutoRefCountObjectPtr da,da2; + if((const ParaMEDMEM::MEDCouplingUMesh *)newMesh) { if((int)ci.size()!=newMesh->getNumberOfNodes()) { da=newMesh->getCellIdsFullyIncludedInNodeIds(&ci[0],&ci[ci.size()]); - mesh2=dynamic_cast(newMesh->buildPartAndReduceNodes(da->getConstPointer(),da->getConstPointer()+da->getNbOfElems(),da2)); + DataArrayInt *tmpp=0; + mesh2=dynamic_cast(newMesh->buildPartAndReduceNodes(da->begin(),da->end(),tmpp)); da2=tmpp; } } else @@ -1917,7 +1927,8 @@ ParaMEDMEM::MEDCouplingFieldDouble *MEDLoaderNS::readFieldDoubleLev2(const char if((int)ci.size()!=mesh->getNumberOfNodes()) { da=mesh->getCellIdsFullyIncludedInNodeIds(&ci[0],&ci[ci.size()]); - mesh2=dynamic_cast(mesh->buildPartAndReduceNodes(da->getConstPointer(),da->getConstPointer()+da->getNbOfElems(),da2)); + DataArrayInt *tmpp=0; + mesh2=dynamic_cast(mesh->buildPartAndReduceNodes(da->begin(),da->end(),tmpp)); da2=tmpp; // int nnodes=mesh2->getNumberOfNodes(); MEDCouplingAutoRefCountObjectPtr da3=DataArrayInt::New(); @@ -1942,32 +1953,47 @@ ParaMEDMEM::MEDCouplingFieldDouble *MEDLoaderNS::readFieldDoubleLev2(const char mesh2->renumberNodes(da2->getConstPointer(),(int)ci.size()); } } - if(da) - da->decrRef(); - if(da2) - da2->decrRef(); + newMesh=mesh2; } - if(newMesh) - newMesh->decrRef(); - newMesh=mesh2; + } + } + else + { + newMesh=const_cast(static_cast(mesh)); mesh->incrRef(); + std::vector types; + for(std::list::const_iterator iter=fieldPerCellType.begin();iter!=fieldPerCellType.end();iter++) + if(std::find(types.begin(),types.end(),(*iter).getType())==types.end()) + types.push_back((*iter).getType()); + for(std::vector::const_iterator it=types.begin();it!=types.end();it++) + { + std::vector cids; + for(std::list::const_iterator iter=fieldPerCellType.begin();iter!=fieldPerCellType.end();iter++) + { + if((*iter).getType()==*it) + { + const std::vector& cellIds=(*iter).getCellIdPerType(); + if(!cellIds.empty()) + std::transform(cellIds.begin(),cellIds.end(),std::back_insert_iterator< std::vector >(cids),std::bind2nd(std::plus(),-1)); + } + } + if(!cids.empty()) + newMesh=newMesh->keepSpecifiedCells(*it,&cids[0],&cids[0]+cids.size()); } } // ParaMEDMEM::MEDCouplingFieldDouble *ret=ParaMEDMEM::MEDCouplingFieldDouble::New(typeOfOutField,ONE_TIME); ret->setName(fieldName); ret->setTime(time,iteration,order); - if(newMesh) + ParaMEDMEM::DataArrayDouble *arr=buildArrayFromRawData(fieldPerCellType,infos); + ret->setArray(arr); + arr->decrRef(); + if((const ParaMEDMEM::MEDCouplingUMesh *)newMesh) { newMesh->setName(mName.c_str());//retrieving mesh name to avoid renaming due to mesh restriction in case of profile. ret->setMesh(newMesh); - newMesh->decrRef(); } else ret->setMesh(mesh); - ParaMEDMEM::DataArrayDouble *arr=buildArrayFromRawData(fieldPerCellType,infos); - ret->setArray(arr); - arr->decrRef(); - // if(typeOfOutField==ON_GAUSS_PT) fillGaussDataOnField(fileName,fieldPerCellType,ret); if(cellRenum) @@ -2578,10 +2604,10 @@ void MEDLoaderNS::prepareCellFieldDoubleForWriting(const ParaMEDMEM::MEDCoupling curType=(INTERP_KERNEL::NormalizedCellType)conn[*pt]; const int *pt2=std::find_if(pt+1,connI+nbOfCells,ConnReaderML(conn,(int)curType)); if(!cellIdsPerType) - split.push_back(MEDLoader::MEDFieldDoublePerCellType(curType,0,nbComp,pt2-pt,0,0)); + split.push_back(MEDLoader::MEDFieldDoublePerCellType(curType,0,nbComp,pt2-pt,1,0,0)); else { - split.push_back(MEDLoader::MEDFieldDoublePerCellType(curType,0,nbComp,pt2-pt,wCellIdsPT,0)); + split.push_back(MEDLoader::MEDFieldDoublePerCellType(curType,0,nbComp,pt2-pt,1,wCellIdsPT,0)); wCellIdsPT+=std::distance(pt,pt2); } pt=pt2; @@ -2637,7 +2663,7 @@ void MEDLoaderNS::writeFieldTryingToFitExistingMesh(const char *fileName, const throw INTERP_KERNEL::Exception(oss.str().c_str()); } MEDCouplingAutoRefCountObjectPtr m=MEDLoader::ReadUMeshFromFile(fileName,f->getMesh()->getName(),f2); - MEDCouplingAutoRefCountObjectPtr m2=MEDCouplingUMesh::MergeUMeshes(m,(MEDCouplingUMesh *)f->getMesh()); + MEDCouplingAutoRefCountObjectPtr m2=MEDCouplingUMesh::MergeUMeshes(m,static_cast(f->getMesh())); bool areNodesMerged; int newNbOfNodes; MEDCouplingAutoRefCountObjectPtr da=m2->mergeNodes(MEDLoader::_EPS_FOR_NODE_COMP,areNodesMerged,newNbOfNodes); diff --git a/src/MEDLoader/MEDLoader.hxx b/src/MEDLoader/MEDLoader.hxx index 61536db38..7fd9f43b8 100644 --- a/src/MEDLoader/MEDLoader.hxx +++ b/src/MEDLoader/MEDLoader.hxx @@ -66,17 +66,19 @@ class MEDLOADER_EXPORT MEDLoader class MEDFieldDoublePerCellType { public: - MEDFieldDoublePerCellType(INTERP_KERNEL::NormalizedCellType type, double *values, int ncomp, int ntuple, const int *cellIdPerType, const char *locName); + MEDFieldDoublePerCellType(INTERP_KERNEL::NormalizedCellType type, double *values, int ncomp, int nGeoElt, int nbi, const int *cellIdPerType, const char *locName); INTERP_KERNEL::NormalizedCellType getType() const { return _type; } int getNbComp() const { return _ncomp; } - int getNbOfTuple() const { return _ntuple; } - int getNbOfValues() const { return _ncomp*_ntuple; } + int getNbOfGeoElt() const { return _ngeo_elt; } + int getNbOfTuple() const { return _nbi*_ngeo_elt; } + int getNbOfValues() const { return _ncomp*_nbi*_ngeo_elt; } double *getArray() const { return _values; } const std::string& getLocName() const { return _loc_name; } const std::vector& getCellIdPerType() const { return _cell_id_per_type; } void releaseArray(); private: - int _ntuple; + int _ngeo_elt; + int _nbi; int _ncomp; double *_values; std::string _loc_name; diff --git a/src/MEDLoader/Swig/MEDLoaderTest3.py b/src/MEDLoader/Swig/MEDLoaderTest3.py index 1622b2980..1c7276bec 100644 --- a/src/MEDLoader/Swig/MEDLoaderTest3.py +++ b/src/MEDLoader/Swig/MEDLoaderTest3.py @@ -1432,6 +1432,141 @@ class MEDLoaderTest(unittest.TestCase): # mfd.write(fname,2) pass + + def testGaussWriteOnPfl1(self): + fname="Pyfile49.med" + coords=DataArrayDouble([0.,0.,0.,1.,1.,1.,1.,0.,0.,0.5,0.5,1.,1.,0.5,0.5,0.],8,2) + mQ8=MEDCouplingUMesh("",2) ; mQ8.setCoords(coords) + mQ8.allocateCells(1) + mQ8.insertNextCell(NORM_QUAD8,range(8)) + mQ8.finishInsertingCells() + mQ4=MEDCouplingUMesh("",2) ; mQ4.setCoords(coords) + mQ4.allocateCells(1) + mQ4.insertNextCell(NORM_QUAD4,range(4)) + mQ4.finishInsertingCells() + mT3=MEDCouplingUMesh("",2) ; mT3.setCoords(coords) + mT3.allocateCells(1) + mT3.insertNextCell(NORM_TRI3,range(3)) + mT3.finishInsertingCells() + + tr=[[0.,4.],[2.,4.],[4.,4.],[6.,4.],[8.,4.],[10.,4.],[12.,4.],[14.,4.],[16.,4.],[18.,4.],[20.,4.],[0.,0.],[2.,0.], [0.,2.],[2.,2.],[4.,2.],[6.,2.],[8.,2.],[10.,2.],[12.,2.]] + ms=11*[mT3]+2*[mQ4]+7*[mQ8] + ms[:]=(elt.deepCpy() for elt in ms) + for m,t in zip(ms,tr): + d=m.getCoords() ; d+= t + pass + m=MEDCouplingUMesh.MergeUMeshes(ms) + m.setName("mesh") + m2=m[:13] ; m2.setName(m.getName()) + ### Use case 1 : Pfl on all tri3 and on all quad4. If we were on CELLS or GAUSS_NE no pfl were needed. But here 2 discs in tri3. + ### So here 2 pfls will be created (pfl_TRI3_loc_0 and pfl_TRI3_loc_1) + f=MEDCouplingFieldDouble.New(ON_GAUSS_PT,ONE_TIME) + f.setMesh(m2) + f.setTime(4.5,1,2) + da=DataArrayDouble(34) ; da.iota(3.) + f.setArray(da) + f.setName("fieldCellOnPflWithoutPfl") + f.setGaussLocalizationOnCells([0,1,2,3,4,5,6,7,8],[0.,0.,1.,0.,1.,1.],[0.3,0.3,0.7,0.7],[0.8,0.2]) + f.setGaussLocalizationOnCells([9,10],[0.,0.,1.,0.,1.,1.],[0.3,0.3,0.7,0.7,0.8,0.8],[0.8,0.07,0.13]) + f.setGaussLocalizationOnCells([11,12],[0.,0.,1.,0.,1.,1.,0.,1.],[0.3,0.3,0.7,0.7,0.8,0.8,0.8,0.8,0.8,0.8],[0.8,0.07,0.1,0.01,0.02]) + f.checkCoherency() + # + mm=MEDFileUMesh() + mm.setMeshAtLevel(0,m) + mm.write(fname,2) + # + f1ts=MEDFileField1TS.New() + pfl=DataArrayInt(range(13)) ; pfl.setName("pfl") + f1ts.setFieldProfile(f,mm,0,pfl) + f1ts.write(fname,0) + # + self.assertEqual(f1ts.getPfls(),('pfl_NORM_TRI3_loc_0', 'pfl_NORM_TRI3_loc_1')) + self.assertEqual(f1ts.getPflsReallyUsed(),('pfl_NORM_TRI3_loc_0', 'pfl_NORM_TRI3_loc_1')) + da1=DataArrayInt([0,1,2,3,4,5,6,7,8]) ; da1.setName("pfl_NORM_TRI3_loc_0") + self.assertTrue(f1ts.getProfile("pfl_NORM_TRI3_loc_0").isEqual(da1)) + da1=DataArrayInt([9,10]) ; da1.setName("pfl_NORM_TRI3_loc_1") + self.assertTrue(f1ts.getProfile("pfl_NORM_TRI3_loc_1").isEqual(da1)) + self.assertEqual(f1ts.getLocs(),('Loc_fieldCellOnPflWithoutPfl_NORM_TRI3_0', 'Loc_fieldCellOnPflWithoutPfl_NORM_TRI3_1', 'Loc_fieldCellOnPflWithoutPfl_NORM_QUAD4_2')) + self.assertEqual(f1ts.getLocsReallyUsed(),('Loc_fieldCellOnPflWithoutPfl_NORM_TRI3_0', 'Loc_fieldCellOnPflWithoutPfl_NORM_TRI3_1', 'Loc_fieldCellOnPflWithoutPfl_NORM_QUAD4_2')) + # + dataRead=MEDFileData.New(fname) + mRead=dataRead.getMeshes()[0] + f1tsRead=dataRead.getFields()[0][0] + f1tsRead.getFieldOnMeshAtLevel(ON_GAUSS_PT,0,mRead) + f2=f1tsRead.getFieldOnMeshAtLevel(ON_GAUSS_PT,0,mRead) + self.assertTrue(f.isEqual(f2,1e-12,1e-12)) + f2_bis=MEDLoader.ReadFieldGauss(fname,m.getName(),0,f.getName(),f.getTime()[1],f.getTime()[2]) + f2_bis.checkCoherency() + self.assertTrue(f.isEqual(f2_bis,1e-12,1e-12)) + ## Use case 2 : Pfl on part tri3 with 2 disc and on part quad8 with 1 disc + f=MEDCouplingFieldDouble.New(ON_GAUSS_PT,ONE_TIME) + pfl=DataArrayInt([1,2,5,6,8,9,15,16,17,18]) ; pfl.setName("pfl2") + m2=m[pfl] ; m2.setName(m.getName()) + f.setMesh(m2) + f.setTime(4.5,1,2) + da=DataArrayDouble(35) ; da.iota(3.) + f.setArray(da) + f.setName("fieldCellOnPflWithoutPfl2") + f.setGaussLocalizationOnCells([0,1,3],[0.,0.,1.,0.,1.,1.],[0.3,0.3,0.7,0.7],[0.8,0.2]) + f.setGaussLocalizationOnCells([2,4,5],[0.,0.,1.,0.,1.,1.],[0.3,0.3,0.7,0.7,0.8,0.8],[0.8,0.07,0.13]) + f.setGaussLocalizationOnCells([6,7,8,9],[0.,0.,1.,0.,1.,1.,0.,1.,0.5,0.,1.,0.5,0.5,1.,0.,0.5],[0.3,0.3,0.7,0.7,0.8,0.8,0.8,0.8,0.8,0.8],[0.8,0.07,0.1,0.01,0.02]) + f.checkCoherency() + # + mm=MEDFileUMesh() + mm.setMeshAtLevel(0,m) + mm.write(fname,2) + f1ts=MEDFileField1TS.New() + f1ts.setFieldProfile(f,mm,0,pfl) + self.assertEqual(f1ts.getPfls(),('pfl2_NORM_TRI3_loc_0','pfl2_NORM_TRI3_loc_1','pfl2_NORM_QUAD8_loc_2')) + self.assertEqual(f1ts.getProfile("pfl2_NORM_TRI3_loc_0").getValues(),[1,2,6]) + self.assertEqual(f1ts.getProfile("pfl2_NORM_TRI3_loc_1").getValues(),[5,8,9]) + self.assertEqual(f1ts.getProfile("pfl2_NORM_QUAD8_loc_2").getValues(),[2,3,4,5]) + f1ts.write(fname,0) + dataRead=MEDFileData.New(fname) + mRead=dataRead.getMeshes()[0] + f1tsRead=dataRead.getFields()[0][0] + f1tsRead.getFieldOnMeshAtLevel(ON_GAUSS_PT,0,mRead) + f3=f1tsRead.getFieldOnMeshAtLevel(ON_GAUSS_PT,0,mRead) + f3.renumberCells([0,1,3,2,4,5,6,7,8,9]) + self.assertTrue(f.isEqual(f3,1e-12,1e-12)) + f3_bis=MEDLoader.ReadFieldGauss(fname,m.getName(),0,f.getName(),f.getTime()[1],f.getTime()[2]) + f3_bis.renumberCells([0,1,3,2,4,5,6,7,8,9]) + self.assertTrue(f.isEqual(f3_bis,1e-12,1e-12)) + ## Use case 3 : no pfl but creation of pfls due to gauss pts + f=MEDCouplingFieldDouble.New(ON_GAUSS_PT,ONE_TIME) + f.setMesh(m) + f.setTime(4.5,1,2) + da=DataArrayDouble(60) ; da.iota(3.) + f.setArray(da) + f.setName("fieldCellWithoutPfl") + f.setGaussLocalizationOnCells([0,1,2,3,4,5,6,7,8],[0.,0.,1.,0.,1.,1.],[0.3,0.3,0.7,0.7],[0.8,0.2]) + f.setGaussLocalizationOnCells([9,10],[0.,0.,1.,0.,1.,1.],[0.3,0.3,0.7,0.7,0.8,0.8],[0.8,0.07,0.13]) + f.setGaussLocalizationOnCells([11,12],[0.,0.,1.,0.,1.,1.,0.,1.],[0.3,0.3,0.7,0.7,0.8,0.8,0.8,0.8,0.8,0.8],[0.8,0.07,0.1,0.01,0.02]) + f.setGaussLocalizationOnCells([13,14,15,17,18],[0.,0.,1.,0.,1.,1.,0.,1.,0.5,0.,1.,0.5,0.5,1.,0.,0.5],[0.3,0.3,0.7,0.7,0.8,0.8,0.8,0.8],[0.8,0.1,0.03,0.07]) + f.setGaussLocalizationOnCells([16,19],[0.,0.,1.,0.,1.,1.,0.,1.,0.5,0.,1.,0.5,0.5,1.,0.,0.5],[0.3,0.3,0.7,0.7,0.8,0.8],[0.8,0.1,0.1]) + f.checkCoherency() + mm=MEDFileUMesh() + mm.setMeshAtLevel(0,m) + f1ts=MEDFileField1TS.New() + f1ts.setFieldNoProfileSBT(f) + self.assertEqual(f1ts.getPfls(),('Pfl_fieldCellWithoutPfl_NORM_TRI3_0','Pfl_fieldCellWithoutPfl_NORM_TRI3_1','Pfl_fieldCellWithoutPfl_NORM_QUAD8_3','Pfl_fieldCellWithoutPfl_NORM_QUAD8_4')) + self.assertEqual(f1ts.getProfile("Pfl_fieldCellWithoutPfl_NORM_TRI3_0").getValues(),[0,1,2,3,4,5,6,7,8]) + self.assertEqual(f1ts.getProfile("Pfl_fieldCellWithoutPfl_NORM_TRI3_1").getValues(),[9,10]) + self.assertEqual(f1ts.getProfile("Pfl_fieldCellWithoutPfl_NORM_QUAD8_3").getValues(),[0,1,2,4,5]) + self.assertEqual(f1ts.getProfile("Pfl_fieldCellWithoutPfl_NORM_QUAD8_4").getValues(),[3,6]) + mm.write(fname,2) + f1ts.write(fname,0) + # + dataRead=MEDFileData.New(fname) + mRead=dataRead.getMeshes()[0] + f1tsRead=dataRead.getFields()[0][0] + f3=f1tsRead.getFieldOnMeshAtLevel(ON_GAUSS_PT,0,mRead) + f3.renumberCells([0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,17,18,16,19]) + self.assertTrue(f.isEqual(f3,1e-12,1e-12)) + f3_bis=MEDLoader.ReadFieldGauss(fname,m.getName(),0,f.getName(),f.getTime()[1],f.getTime()[2]) + f3_bis.renumberCells([0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,17,18,16,19]) + self.assertTrue(f.isEqual(f3_bis,1e-12,1e-12)) + pass pass unittest.main()