From ddab81eebca8daf4d4abe656b6eae82c08f5cb7e Mon Sep 17 00:00:00 2001 From: ageay Date: Fri, 30 Sep 2011 09:28:09 +0000 Subject: [PATCH] Introduction of 3D/2D interpolation into MEDCouplingRemapper. --- src/MEDCoupling/MEDCouplingRemapper.cxx | 62 +++++++++++++++++++ .../Test/MEDCouplingRemapperTest.cxx | 37 +++++++++++ 2 files changed, 99 insertions(+) diff --git a/src/MEDCoupling/MEDCouplingRemapper.cxx b/src/MEDCoupling/MEDCouplingRemapper.cxx index 8c9f3ee40..9e675e370 100644 --- a/src/MEDCoupling/MEDCouplingRemapper.cxx +++ b/src/MEDCoupling/MEDCouplingRemapper.cxx @@ -31,6 +31,7 @@ #include "Interpolation3D.txx" #include "Interpolation3DSurf.hxx" #include "Interpolation2D1D.txx" +#include "Interpolation3D2D.txx" using namespace ParaMEDMEM; @@ -251,6 +252,17 @@ int MEDCouplingRemapper::prepareUU(const char *method) throw(INTERP_KERNEL::Exce nbCols=interpolation.interpolateMeshes(target_mesh_wrapper,source_mesh_wrapper,matrixTmp,method); reverseMatrix(matrixTmp,nbCols,_matrix); nbCols=matrixTmp.size(); + INTERP_KERNEL::Interpolation2D1D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces(); + if(!duplicateFaces.empty()) + { + std::ostringstream oss; oss << "An unexpected situation happend ! For the following 1D Cells are part of edges shared by 2D cells :\n"; + for(std::map >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++) + { + oss << "1D Cell #" << (*it).first << " is part of common edge of following 2D cells ids : "; + std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator(oss," ")); + oss << std::endl; + } + } } } else if(srcMeshDim==1 && trgMeshDim==2 && srcSpaceDim==2) @@ -271,6 +283,56 @@ int MEDCouplingRemapper::prepareUU(const char *method) throw(INTERP_KERNEL::Exce MEDCouplingNormalizedUnstructuredMesh<2,2> target_mesh_wrapper(target_mesh); INTERP_KERNEL::Interpolation2D1D interpolation(*this); nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method); + INTERP_KERNEL::Interpolation2D1D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces(); + if(!duplicateFaces.empty()) + { + std::ostringstream oss; oss << "An unexpected situation happend ! For the following 1D Cells are part of edges shared by 2D cells :\n"; + for(std::map >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++) + { + oss << "1D Cell #" << (*it).first << " is part of common edge of following 2D cells ids : "; + std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator(oss," ")); + oss << std::endl; + } + } + } + } + else if(srcMeshDim==2 && trgMeshDim==3 && srcSpaceDim==3) + { + MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh); + MEDCouplingNormalizedUnstructuredMesh<3,3> target_mesh_wrapper(target_mesh); + INTERP_KERNEL::Interpolation3D2D interpolation(*this); + nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method); + INTERP_KERNEL::Interpolation3D2D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces(); + if(!duplicateFaces.empty()) + { + std::ostringstream oss; oss << "An unexpected situation happend ! For the following 2D Cells are part of edges shared by 3D cells :\n"; + for(std::map >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++) + { + oss << "2D Cell #" << (*it).first << " is part of common face of following 3D cells ids : "; + std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator(oss," ")); + oss << std::endl; + } + } + } + else if(srcMeshDim==3 && trgMeshDim==2 && srcSpaceDim==3) + { + MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh); + MEDCouplingNormalizedUnstructuredMesh<3,3> target_mesh_wrapper(target_mesh); + INTERP_KERNEL::Interpolation3D2D interpolation(*this); + std::vector > matrixTmp; + nbCols=interpolation.interpolateMeshes(target_mesh_wrapper,source_mesh_wrapper,matrixTmp,method); + reverseMatrix(matrixTmp,nbCols,_matrix); + nbCols=matrixTmp.size(); + INTERP_KERNEL::Interpolation3D2D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces(); + if(!duplicateFaces.empty()) + { + std::ostringstream oss; oss << "An unexpected situation happend ! For the following 2D Cells are part of edges shared by 3D cells :\n"; + for(std::map >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++) + { + oss << "2D Cell #" << (*it).first << " is part of common face of following 3D cells ids : "; + std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator(oss," ")); + oss << std::endl; + } } } else if(trgMeshDim==-1) diff --git a/src/MEDCoupling/Test/MEDCouplingRemapperTest.cxx b/src/MEDCoupling/Test/MEDCouplingRemapperTest.cxx index 58d2de4f6..3bbb7e2e4 100644 --- a/src/MEDCoupling/Test/MEDCouplingRemapperTest.cxx +++ b/src/MEDCoupling/Test/MEDCouplingRemapperTest.cxx @@ -687,6 +687,43 @@ void MEDCouplingRemapperTest::testMultiDimCombi() trgField->decrRef(); sourceMesh->decrRef(); targetMesh->decrRef(); + //------------- 2D -> 3D + sourceMesh=MEDCouplingBasicsTest::build3D2DSourceMesh(); + targetMesh=MEDCouplingBasicsTest::build3D2DTargetMesh(); + remapper.setIntersectionType(INTERP_KERNEL::Triangulation); + CPPUNIT_ASSERT_EQUAL(1,remapper.prepare(sourceMesh,targetMesh,"P0P0")); + srcField=MEDCouplingFieldDouble::New(ON_CELLS); + srcField->setNature(ConservativeVolumic); + srcField->setMesh(sourceMesh); + array=DataArrayDouble::New(); + array->alloc(7,1); array->iota(2.); + srcField->setArray(array); array->decrRef(); + trgField=remapper.transferField(srcField,4.57); + const double valuesExpected12[3]={5.70909090909091, 6.08362715128042, 6.92857142857143}; + CPPUNIT_ASSERT_EQUAL(3,trgField->getNumberOfTuples()); + CPPUNIT_ASSERT_EQUAL(1,trgField->getNumberOfComponents()); + for(int i=0;i<3;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected12[i],trgField->getIJ(i,0),1e-13); + srcField->decrRef(); + trgField->decrRef(); + //------------- 3D -> 2D + CPPUNIT_ASSERT_EQUAL(1,remapper.prepare(targetMesh,sourceMesh,"P0P0")); + srcField=MEDCouplingFieldDouble::New(ON_CELLS); + srcField->setNature(ConservativeVolumic); + srcField->setMesh(targetMesh); + array=DataArrayDouble::New(); + array->alloc(3,1); array->iota(2.); + srcField->setArray(array); array->decrRef(); + trgField=remapper.transferField(srcField,4.57); + const double valuesExpected13[7]={3., 4., 2.5, 2.909090909090909, 2., 3.5, 3.3571428571428572}; + CPPUNIT_ASSERT_EQUAL(7,trgField->getNumberOfTuples()); + CPPUNIT_ASSERT_EQUAL(1,trgField->getNumberOfComponents()); + for(int i=0;i<7;i++) + CPPUNIT_ASSERT_DOUBLES_EQUAL(valuesExpected13[i],trgField->getIJ(i,0),1e-13); + srcField->decrRef(); + trgField->decrRef(); + sourceMesh->decrRef(); + targetMesh->decrRef(); } void MEDCouplingRemapperTest::testNatureOfField() -- 2.39.2