From e0a6d839ef5bc27fc74ba5c02b9e43b745023eaf Mon Sep 17 00:00:00 2001 From: Anthony Geay Date: Mon, 22 Aug 2022 11:40:19 +0200 Subject: [PATCH] WIP --- .../reference/fields/discretization.dox | 1 + src/INTERP_KERNEL/BBTreeStandAlone.txx | 1 + src/INTERP_KERNEL/InterpolationHelper.txx | 1 + src/MEDCoupling/CMakeLists.txt | 1 + .../MEDCouplingFieldDiscretization.cxx | 57 ++++----- .../MEDCouplingFieldDiscretization.hxx | 86 +++++++------ ...EDCouplingFieldDiscretizationOnNodesFE.cxx | 118 ++++++++++++++++++ ...EDCouplingFieldDiscretizationOnNodesFE.hxx | 49 ++++++++ src/MEDCoupling/MEDCouplingFieldT.txx | 6 +- .../MEDCouplingNormalizedUnstructuredMesh.hxx | 2 +- src/MEDCoupling/MEDCouplingRefCountObject.hxx | 3 +- src/MEDCoupling_Swig/MEDCouplingCommon.i | 3 +- 12 files changed, 252 insertions(+), 76 deletions(-) create mode 100755 src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.cxx create mode 100644 src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.hxx diff --git a/doc/developer/doxygen/doxfiles/reference/fields/discretization.dox b/doc/developer/doxygen/doxfiles/reference/fields/discretization.dox index abb5dd13b..31f4d4b32 100644 --- a/doc/developer/doxygen/doxfiles/reference/fields/discretization.dox +++ b/doc/developer/doxygen/doxfiles/reference/fields/discretization.dox @@ -21,6 +21,7 @@ A field can be supported by: - Gauss points: built with \ref MEDCoupling::TypeOfField "ON_GAUSS_PT" - Gauss points on nodes per element: built with \ref MEDCoupling::TypeOfField "ON_GAUSS_NE" - Kriging points: built with \ref MEDCoupling::TypeOfField "ON_NODES_KR" + - On points using standard FE description : built with \ref MEDCoupling::TypeOfField "ON_NODES_FE" The spatial discretization is at the center of the \ref interpolation "interpolation" mechanisms, since one of the main interpolation parameter is indeed specifying from which source discretization diff --git a/src/INTERP_KERNEL/BBTreeStandAlone.txx b/src/INTERP_KERNEL/BBTreeStandAlone.txx index 73b28580c..00ea437df 100644 --- a/src/INTERP_KERNEL/BBTreeStandAlone.txx +++ b/src/INTERP_KERNEL/BBTreeStandAlone.txx @@ -34,4 +34,5 @@ private: public: BBTreeStandAlone(std::unique_ptr&& bbs, ConnType nbelems, double epsilon=BBTREE_DFT_EPSILON):_bbox(std::move(bbs)),_effective(_bbox.get(),nullptr,0,nbelems,epsilon) { } void getIntersectingElems(const double* bb, std::vector& elems) const { _effective.getIntersectingElems(bb,elems); } + void getElementsAroundPoint(const double* xx, std::vector& elems) const { _effective.getElementsAroundPoint(xx,elems); } }; diff --git a/src/INTERP_KERNEL/InterpolationHelper.txx b/src/INTERP_KERNEL/InterpolationHelper.txx index 4fcde3237..48ca55d1f 100755 --- a/src/INTERP_KERNEL/InterpolationHelper.txx +++ b/src/INTERP_KERNEL/InterpolationHelper.txx @@ -21,6 +21,7 @@ #pragma once #include "BBTreeStandAlone.txx" +#include "MeshElement.txx" #include "Log.hxx" #include diff --git a/src/MEDCoupling/CMakeLists.txt b/src/MEDCoupling/CMakeLists.txt index 8a18a7a55..8dd2592a6 100644 --- a/src/MEDCoupling/CMakeLists.txt +++ b/src/MEDCoupling/CMakeLists.txt @@ -61,6 +61,7 @@ SET(medcoupling_SOURCES MEDCouplingStructuredMesh.cxx MEDCouplingTimeDiscretization.cxx MEDCouplingFieldDiscretization.cxx + MEDCouplingFieldDiscretizationOnNodesFE.cxx MEDCouplingRefCountObject.cxx MEDCouplingPointSet.cxx MEDCouplingFieldTemplate.cxx diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx index 693f24b30..d7cf86b66 100755 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.cxx @@ -19,6 +19,7 @@ // Author : Anthony Geay (EDF R&D) #include "MEDCouplingFieldDiscretization.hxx" +#include "MEDCouplingFieldDiscretizationOnNodesFE.hxx" #include "MEDCouplingCMesh.hxx" #include "MEDCouplingUMesh.hxx" #include "MEDCouplingFieldDouble.hxx" @@ -44,26 +45,16 @@ const double MEDCouplingFieldDiscretization::DFLT_PRECISION=1.e-12; const char MEDCouplingFieldDiscretizationP0::REPR[]="P0"; -const TypeOfField MEDCouplingFieldDiscretizationP0::TYPE=ON_CELLS; - const char MEDCouplingFieldDiscretizationP1::REPR[]="P1"; -const TypeOfField MEDCouplingFieldDiscretizationP1::TYPE=ON_NODES; - const mcIdType MEDCouplingFieldDiscretizationPerCell::DFT_INVALID_LOCID_VALUE=-1; const char MEDCouplingFieldDiscretizationGauss::REPR[]="GAUSS"; -const TypeOfField MEDCouplingFieldDiscretizationGauss::TYPE=ON_GAUSS_PT; - const char MEDCouplingFieldDiscretizationGaussNE::REPR[]="GSSNE"; -const TypeOfField MEDCouplingFieldDiscretizationGaussNE::TYPE=ON_GAUSS_NE; - const char MEDCouplingFieldDiscretizationKriging::REPR[]="KRIGING"; -const TypeOfField MEDCouplingFieldDiscretizationKriging::TYPE=ON_NODES_KR; - // doc is here http://www.code-aster.org/V2/doc/default/fr/man_r/r3/r3.01.01.pdf const double MEDCouplingFieldDiscretizationGaussNE::FGP_POINT1[1]={0.}; const double MEDCouplingFieldDiscretizationGaussNE::FGP_SEG2[2]={1.,1.}; @@ -142,6 +133,8 @@ MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretization::New(TypeOfField return new MEDCouplingFieldDiscretizationGaussNE; case MEDCouplingFieldDiscretizationKriging::TYPE: return new MEDCouplingFieldDiscretizationKriging; + case MEDCouplingFieldDiscretizationOnNodesFE::TYPE: + return new MEDCouplingFieldDiscretizationOnNodesFE; default: throw INTERP_KERNEL::Exception("Chosen discretization is not implemented yet."); } @@ -159,22 +152,30 @@ TypeOfField MEDCouplingFieldDiscretization::GetTypeOfFieldFromStringRepr(const s return MEDCouplingFieldDiscretizationGaussNE::TYPE; if(repr==MEDCouplingFieldDiscretizationKriging::REPR) return MEDCouplingFieldDiscretizationKriging::TYPE; + if(repr==MEDCouplingFieldDiscretizationOnNodesFE::REPR) + return MEDCouplingFieldDiscretizationOnNodesFE::TYPE; throw INTERP_KERNEL::Exception("Representation does not match with any field discretization !"); } std::string MEDCouplingFieldDiscretization::GetTypeOfFieldRepr(TypeOfField type) { - if(type==MEDCouplingFieldDiscretizationP0::TYPE) - return MEDCouplingFieldDiscretizationP0::REPR; - if(type==MEDCouplingFieldDiscretizationP1::TYPE) - return MEDCouplingFieldDiscretizationP1::REPR; - if(type==MEDCouplingFieldDiscretizationGauss::TYPE) - return MEDCouplingFieldDiscretizationGauss::REPR; - if(type==MEDCouplingFieldDiscretizationGaussNE::TYPE) - return MEDCouplingFieldDiscretizationGaussNE::REPR; - if(type==MEDCouplingFieldDiscretizationKriging::TYPE) - return MEDCouplingFieldDiscretizationKriging::REPR; - throw INTERP_KERNEL::Exception("GetTypeOfFieldRepr : Representation does not match with any field discretization !"); + switch(type) + { + case MEDCouplingFieldDiscretizationP0::TYPE: + return MEDCouplingFieldDiscretizationP0::REPR; + case MEDCouplingFieldDiscretizationP1::TYPE: + return MEDCouplingFieldDiscretizationP1::REPR; + case MEDCouplingFieldDiscretizationGauss::TYPE: + return MEDCouplingFieldDiscretizationGauss::REPR; + case MEDCouplingFieldDiscretizationGaussNE::TYPE: + return MEDCouplingFieldDiscretizationGaussNE::REPR; + case MEDCouplingFieldDiscretizationKriging::TYPE: + return MEDCouplingFieldDiscretizationKriging::REPR; + case MEDCouplingFieldDiscretizationOnNodesFE::TYPE: + return MEDCouplingFieldDiscretizationOnNodesFE::REPR; + default: + throw INTERP_KERNEL::Exception("GetTypeOfFieldRepr : Representation does not match with any field discretization !"); + } } bool MEDCouplingFieldDiscretization::isEqual(const MEDCouplingFieldDiscretization *other, double eps) const @@ -480,20 +481,6 @@ void MEDCouplingFieldDiscretization::RenumberEntitiesFromN2OArr(const mcIdType * } } -template -MCAuto MEDCouplingFieldDiscretization::EasyAggregate(std::vector& fds) -{ - if(fds.empty()) - throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretization::aggregate : input array is empty"); - for(const MEDCouplingFieldDiscretization * it : fds) - { - const FIELD_DISC *itc(dynamic_cast(it)); - if(!itc) - throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretization::aggregate : same field discretization expected for all input discretizations !"); - } - return fds[0]->clone(); -} - MEDCouplingFieldDiscretization::~MEDCouplingFieldDiscretization() { } diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx index e1004c571..2eab0a945 100644 --- a/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx +++ b/src/MEDCoupling/MEDCouplingFieldDiscretization.hxx @@ -124,10 +124,10 @@ namespace MEDCoupling public: MEDCOUPLING_EXPORT TypeOfField getEnum() const; MEDCOUPLING_EXPORT std::string getClassName() const override { return std::string("MEDCouplingFieldDiscretizationP0"); } - MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const; + MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const override; MEDCOUPLING_EXPORT std::string getStringRepr() const; MEDCOUPLING_EXPORT const char *getRepr() const; - MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const; + MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const override; MEDCOUPLING_EXPORT mcIdType getNumberOfTuplesExpectedRegardingCode(const std::vector& code, const std::vector& idsPerType) const; MEDCOUPLING_EXPORT mcIdType getNumberOfTuples(const MEDCouplingMesh *mesh) const; MEDCOUPLING_EXPORT mcIdType getNumberOfMeshPlaces(const MEDCouplingMesh *mesh) const; @@ -135,14 +135,14 @@ namespace MEDCoupling MEDCOUPLING_EXPORT void renumberArraysForCell(const MEDCouplingMesh *mesh, const std::vector& arrays, const mcIdType *old2NewBg, bool check); MEDCOUPLING_EXPORT DataArrayDouble *getLocalizationOfDiscValues(const MEDCouplingMesh *mesh) const; - MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const; + MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const override; MEDCOUPLING_EXPORT void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const mcIdType *tupleIdsBg, const mcIdType *tupleIdsEnd, DataArrayIdType *&cellRestriction, DataArrayIdType *&trueTupleRestriction) const; MEDCOUPLING_EXPORT void checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArray *da) const; - MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const; - MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; + MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const override; + MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const override; MEDCOUPLING_EXPORT void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, mcIdType i, mcIdType j, mcIdType k, double *res) const; - MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const; + MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const override; MEDCOUPLING_EXPORT void renumberValuesOnNodes(double epsOnVals, const mcIdType *old2New, mcIdType newNbOfNodes, DataArrayDouble *arr) const; MEDCOUPLING_EXPORT void renumberValuesOnCells(double epsOnVals, const MEDCouplingMesh *mesh, const mcIdType *old2New, mcIdType newSz, DataArrayDouble *arr) const; MEDCOUPLING_EXPORT void renumberValuesOnCellsR(const MEDCouplingMesh *mesh, const mcIdType *new2old, mcIdType newSz, DataArrayDouble *arr) const; @@ -153,7 +153,7 @@ namespace MEDCoupling MEDCOUPLING_EXPORT void reprQuickOverview(std::ostream& stream) const; public: static const char REPR[]; - static const TypeOfField TYPE; + static constexpr TypeOfField TYPE = ON_CELLS; }; class MEDCouplingFieldDiscretizationOnNodes : public MEDCouplingFieldDiscretization @@ -184,19 +184,19 @@ namespace MEDCoupling public: MEDCOUPLING_EXPORT TypeOfField getEnum() const; MEDCOUPLING_EXPORT std::string getClassName() const override { return std::string("MEDCouplingFieldDiscretizationP1"); } - MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const; + MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const override; MEDCOUPLING_EXPORT std::string getStringRepr() const; MEDCOUPLING_EXPORT const char *getRepr() const; - MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const; - MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const; - MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const; - MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; - MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const; + MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const override; + MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const override; + MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const override; + MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const override; + MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const override; MEDCOUPLING_EXPORT void reprQuickOverview(std::ostream& stream) const; MEDCOUPLING_EXPORT MCAuto aggregate(std::vector& fds) const override; public: static const char REPR[]; - static const TypeOfField TYPE; + static constexpr TypeOfField TYPE = ON_NODES; protected: MEDCOUPLING_EXPORT void getValueInCell(const MEDCouplingMesh *mesh, mcIdType cellId, const DataArrayDouble *arr, const double *loc, double *res) const; }; @@ -222,7 +222,7 @@ namespace MEDCoupling std::size_t getHeapMemorySizeWithoutChildren() const; std::vector getDirectChildrenWithNull() const; void checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArray *da) const; - bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const; + bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const override; bool isEqualWithoutConsideringStr(const MEDCouplingFieldDiscretization *other, double eps) const; void renumberCells(const mcIdType *old2NewBg, bool check); protected: @@ -238,9 +238,9 @@ namespace MEDCoupling MEDCOUPLING_EXPORT MEDCouplingFieldDiscretizationGauss(); MEDCOUPLING_EXPORT TypeOfField getEnum() const; MEDCOUPLING_EXPORT std::string getClassName() const override { return std::string("MEDCouplingFieldDiscretizationGauss"); } - MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const; + MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const override; MEDCOUPLING_EXPORT bool isEqualWithoutConsideringStr(const MEDCouplingFieldDiscretization *other, double eps) const; - MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const; + MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const override; MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clonePart(const mcIdType *startCellIds, const mcIdType *endCellIds) const; MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clonePartRange(mcIdType beginCellIds, mcIdType endCellIds, mcIdType stepCellIds) const; MEDCOUPLING_EXPORT std::string getStringRepr() const; @@ -255,7 +255,7 @@ namespace MEDCoupling MEDCOUPLING_EXPORT DataArrayDouble *getLocalizationOfDiscValues(const MEDCouplingMesh *mesh) const; MEDCOUPLING_EXPORT void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const mcIdType *tupleIdsBg, const mcIdType *tupleIdsEnd, DataArrayIdType *&cellRestriction, DataArrayIdType *&trueTupleRestriction) const; - MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const; + MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const override; MEDCOUPLING_EXPORT void getTinySerializationIntInformation(std::vector& tinyInfo) const; MEDCOUPLING_EXPORT void getTinySerializationDbleInformation(std::vector& tinyInfo) const; MEDCOUPLING_EXPORT void finishUnserialization(const std::vector& tinyInfo); @@ -264,10 +264,10 @@ namespace MEDCoupling MEDCOUPLING_EXPORT void checkForUnserialization(const std::vector& tinyInfo, const DataArrayIdType *arr); MEDCOUPLING_EXPORT double getIJK(const MEDCouplingMesh *mesh, const DataArrayDouble *da, mcIdType cellId, mcIdType nodeIdInCell, int compoId) const; MEDCOUPLING_EXPORT void checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArray *da) const; - MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const; - MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; + MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const override; + MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const override; MEDCOUPLING_EXPORT void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, mcIdType i, mcIdType j, mcIdType k, double *res) const; - MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const; + MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const override; MEDCOUPLING_EXPORT MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const mcIdType *start, const mcIdType *end, DataArrayIdType *&di) const; MEDCOUPLING_EXPORT MEDCouplingMesh *buildSubMeshDataRange(const MEDCouplingMesh *mesh, mcIdType beginCellIds, mcIdType endCellIds, mcIdType stepCellIds, mcIdType& beginOut, mcIdType& endOut, mcIdType& stepOut, DataArrayIdType *&di) const; MEDCOUPLING_EXPORT DataArrayIdType *computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const mcIdType *startCellIds, const mcIdType *endCellIds) const; @@ -301,7 +301,7 @@ namespace MEDCoupling void commonUnserialization(const std::vector& tinyInfo); public: static const char REPR[]; - static const TypeOfField TYPE; + static constexpr TypeOfField TYPE = ON_GAUSS_PT; private: std::vector _loc; }; @@ -315,10 +315,10 @@ namespace MEDCoupling MEDCOUPLING_EXPORT MEDCouplingFieldDiscretizationGaussNE(); MEDCOUPLING_EXPORT TypeOfField getEnum() const; MEDCOUPLING_EXPORT std::string getClassName() const override { return std::string("MEDCouplingFieldDiscretizationGaussNE"); } - MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const; + MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const override; MEDCOUPLING_EXPORT std::string getStringRepr() const; MEDCOUPLING_EXPORT const char *getRepr() const; - MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const; + MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const override; MEDCOUPLING_EXPORT mcIdType getNumberOfTuplesExpectedRegardingCode(const std::vector& code, const std::vector& idsPerType) const; MEDCOUPLING_EXPORT mcIdType getNumberOfTuples(const MEDCouplingMesh *mesh) const; MEDCOUPLING_EXPORT mcIdType getNumberOfMeshPlaces(const MEDCouplingMesh *mesh) const; @@ -329,13 +329,13 @@ namespace MEDCoupling MEDCOUPLING_EXPORT void integral(const MEDCouplingMesh *mesh, const DataArrayDouble *arr, bool isWAbs, double *res) const; MEDCOUPLING_EXPORT void computeMeshRestrictionFromTupleIds(const MEDCouplingMesh *mesh, const mcIdType *tupleIdsBg, const mcIdType *tupleIdsEnd, DataArrayIdType *&cellRestriction, DataArrayIdType *&trueTupleRestriction) const; - MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const; + MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const override; MEDCOUPLING_EXPORT double getIJK(const MEDCouplingMesh *mesh, const DataArrayDouble *da, mcIdType cellId, mcIdType nodeIdInCell, int compoId) const; MEDCOUPLING_EXPORT void checkCoherencyBetween(const MEDCouplingMesh *mesh, const DataArray *da) const; - MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const; - MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; + MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const override; + MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const override; MEDCOUPLING_EXPORT void getValueOnPos(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, mcIdType i, mcIdType j, mcIdType k, double *res) const; - MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const; + MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const override; MEDCOUPLING_EXPORT MEDCouplingMesh *buildSubMeshData(const MEDCouplingMesh *mesh, const mcIdType *start, const mcIdType *end, DataArrayIdType *&di) const; MEDCOUPLING_EXPORT MEDCouplingMesh *buildSubMeshDataRange(const MEDCouplingMesh *mesh, mcIdType beginCellIds, mcIdType endCellIds, mcIdType stepCellIds, mcIdType& beginOut, mcIdType& endOut, mcIdType& stepOut, DataArrayIdType *&di) const; MEDCOUPLING_EXPORT DataArrayIdType *computeTupleIdsToSelectFromCellIds(const MEDCouplingMesh *mesh, const mcIdType *startCellIds, const mcIdType *endCellIds) const; @@ -351,7 +351,7 @@ namespace MEDCoupling MEDCOUPLING_EXPORT MEDCouplingFieldDiscretizationGaussNE(const MEDCouplingFieldDiscretizationGaussNE& other); public: static const char REPR[]; - static const TypeOfField TYPE; + static constexpr TypeOfField TYPE = ON_GAUSS_NE; static const double FGP_POINT1[1]; static const double FGP_SEG2[2]; static const double FGP_SEG3[3]; @@ -418,13 +418,13 @@ namespace MEDCoupling MEDCOUPLING_EXPORT TypeOfField getEnum() const; MEDCOUPLING_EXPORT std::string getClassName() const override { return std::string("MEDCouplingFieldDiscretizationKriging"); } MEDCOUPLING_EXPORT const char *getRepr() const; - MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const; + MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const override; MEDCOUPLING_EXPORT std::string getStringRepr() const; - MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const; - MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const; - MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const; - MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const; - MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const; + MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const override; + MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const override; + MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const override; + MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const override; + MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const override; MEDCOUPLING_EXPORT void reprQuickOverview(std::ostream& stream) const; MEDCOUPLING_EXPORT MCAuto aggregate(std::vector& fds) const override; public://specific part @@ -440,8 +440,22 @@ namespace MEDCoupling MEDCOUPLING_EXPORT static DataArrayDouble *PerformDriftOfVec(const DataArrayDouble *arr, mcIdType isDrift); public: static const char REPR[]; - static const TypeOfField TYPE; + static constexpr TypeOfField TYPE = ON_NODES_KR; }; + + template + MCAuto MEDCouplingFieldDiscretization::EasyAggregate(std::vector& fds) + { + if(fds.empty()) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretization::aggregate : input array is empty"); + for(const MEDCouplingFieldDiscretization * it : fds) + { + const FIELD_DISC *itc(dynamic_cast(it)); + if(!itc) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretization::aggregate : same field discretization expected for all input discretizations !"); + } + return fds[0]->clone(); + } } #endif diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.cxx b/src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.cxx new file mode 100755 index 000000000..4494bec34 --- /dev/null +++ b/src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.cxx @@ -0,0 +1,118 @@ +// Copyright (C) 2007-2022 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +// Author : Anthony Geay (EDF R&D) + +#include "MEDCouplingFieldDiscretizationOnNodesFE.hxx" +#include "MEDCouplingNormalizedUnstructuredMesh.txx" +#include "MEDCouplingUMesh.hxx" +#include "InterpolationHelper.txx" + +#include + +using namespace MEDCoupling; + +const char MEDCouplingFieldDiscretizationOnNodesFE::REPR[]="FE"; + +std::string MEDCouplingFieldDiscretizationOnNodesFE::getStringRepr() const +{ + return std::string(REPR); +} + +void MEDCouplingFieldDiscretizationOnNodesFE::reprQuickOverview(std::ostream& stream) const +{ + stream << "NodeFE spatial discretization."; +} + +MCAuto MEDCouplingFieldDiscretizationOnNodesFE::aggregate(std::vector& fds) const +{ + return EasyAggregate(fds); +} + +bool MEDCouplingFieldDiscretizationOnNodesFE::isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const +{ + if(!other) + { + reason="other spatial discretization is NULL, and this spatial discretization (Node FE) is defined."; + return false; + } + const MEDCouplingFieldDiscretizationOnNodesFE *otherC=dynamic_cast(other); + bool ret=otherC!=0; + if(!ret) + reason="Spatial discrtization of this is ON_NODES_FE, which is not the case of other."; + return ret; +} + +/*! + * This method is simply called by MEDCouplingFieldDiscretization::deepCopy. It performs the deep copy of \a this. + * + * \sa MEDCouplingFieldDiscretization::deepCopy. + */ +MEDCouplingFieldDiscretization *MEDCouplingFieldDiscretizationOnNodesFE::clone() const +{ + return new MEDCouplingFieldDiscretizationOnNodesFE; +} + +void MEDCouplingFieldDiscretizationOnNodesFE::checkCompatibilityWithNature(NatureOfField nat) const +{ + if(nat!=IntensiveMaximum) + throw INTERP_KERNEL::Exception("Invalid nature for NodeFE field : expected IntensiveMaximum !"); +} + +MEDCouplingFieldDouble *MEDCouplingFieldDiscretizationOnNodesFE::getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const +{ + if(!mesh) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationOnNodesFE::getMeasureField : mesh instance specified is NULL !"); + throw INTERP_KERNEL::Exception("getMeasureField on MEDCouplingFieldDiscretizationOnNodesFE : not implemented yet !"); +} + +void MEDCouplingFieldDiscretizationOnNodesFE::getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const +{ + const MEDCouplingUMesh *umesh( dynamic_cast(mesh) ); + if( !umesh ) + THROW_IK_EXCEPTION("getValueOn : not implemented yet for type != MEDCouplingUMesh !"); + if(umesh->getSpaceDimension() != 3 || umesh->getMeshDimension() != 3) + THROW_IK_EXCEPTION("getValueOn : implemented only for meshes with spacedim == 3 and meshdim == 3 !"); + MEDCouplingNormalizedUnstructuredMesh<3,3> mesh_wrapper(umesh); + BBTreeStandAlone<3,mcIdType> tree( INTERP_KERNEL::BuildBBTree( mesh_wrapper ) ); + std::vector elems; + tree.getElementsAroundPoint(loc,elems); + + /*MCAuto res2(this->getValueOnMulti(arr,mesh,loc,1)); + std::copy(res2->begin(),res2->end(),res);*/ + throw INTERP_KERNEL::Exception("getValueOn on MEDCouplingFieldDiscretizationOnNodesFE : not implemented yet !"); +} + +DataArrayDouble *MEDCouplingFieldDiscretizationOnNodesFE::getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfTargetPoints) const +{ + if(!arr || !arr->isAllocated()) + throw INTERP_KERNEL::Exception("MEDCouplingFieldDiscretizationOnNodesFE::getValueOnMulti : input array is null or not allocated !"); + mcIdType nbOfRows=getNumberOfMeshPlaces(mesh); + if(arr->getNumberOfTuples()!=nbOfRows) + { + THROW_IK_EXCEPTION( "MEDCouplingFieldDiscretizationOnNodesFE::getValueOnMulti : input array does not have correct number of tuples ! Excepted " << nbOfRows << " having " << arr->getNumberOfTuples() << " !") + } + throw INTERP_KERNEL::Exception("getValueOnMulti on MEDCouplingFieldDiscretizationOnNodesFE : not implemented yet !"); + /*mcIdType nbCols(-1); + std::size_t nbCompo=arr->getNumberOfComponents(); + MCAuto m(computeEvaluationMatrixOnGivenPts(mesh,loc,nbOfTargetPoints,nbCols)); + MCAuto ret(DataArrayDouble::New()); + ret->alloc(nbOfTargetPoints,nbCompo); + INTERP_KERNEL::matrixProduct(m->begin(),nbOfTargetPoints,nbCols,arr->begin(),nbOfRows,ToIdType(nbCompo),ret->getPointer()); + return ret.retn();*/ +} diff --git a/src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.hxx b/src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.hxx new file mode 100644 index 000000000..04559381b --- /dev/null +++ b/src/MEDCoupling/MEDCouplingFieldDiscretizationOnNodesFE.hxx @@ -0,0 +1,49 @@ +// Copyright (C) 2022 CEA/DEN, EDF R&D +// +// This library is free software; you can redistribute it and/or +// modify it under the terms of the GNU Lesser General Public +// License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// +// This library is distributed in the hope that it will be useful, +// but WITHOUT ANY WARRANTY; without even the implied warranty of +// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +// Lesser General Public License for more details. +// +// You should have received a copy of the GNU Lesser General Public +// License along with this library; if not, write to the Free Software +// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA +// +// See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com +// +// Author : Anthony Geay (EDF R&D) + +#pragma once + +#include "MEDCouplingFieldDiscretization.hxx" + +namespace MEDCoupling +{ + /*! + * Class in charge to implement FE functions with shape functions + */ + class MEDCouplingFieldDiscretizationOnNodesFE : public MEDCouplingFieldDiscretizationOnNodes + { + public: + MEDCOUPLING_EXPORT TypeOfField getEnum() const override { return TYPE; } + MEDCOUPLING_EXPORT std::string getClassName() const override { return std::string("MEDCouplingFieldDiscretizationOnNodesFE"); } + MEDCOUPLING_EXPORT const char *getRepr() const override { return REPR; } + MEDCOUPLING_EXPORT std::string getStringRepr() const override; + MEDCOUPLING_EXPORT void reprQuickOverview(std::ostream& stream) const override; + MEDCOUPLING_EXPORT MCAuto aggregate(std::vector& fds) const override; + MEDCOUPLING_EXPORT bool isEqualIfNotWhy(const MEDCouplingFieldDiscretization *other, double eps, std::string& reason) const override; + MEDCOUPLING_EXPORT MEDCouplingFieldDiscretization *clone() const override; + MEDCOUPLING_EXPORT void checkCompatibilityWithNature(NatureOfField nat) const override; + MEDCOUPLING_EXPORT MEDCouplingFieldDouble *getMeasureField(const MEDCouplingMesh *mesh, bool isAbs) const override; + MEDCOUPLING_EXPORT void getValueOn(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, double *res) const override; + MEDCOUPLING_EXPORT DataArrayDouble *getValueOnMulti(const DataArrayDouble *arr, const MEDCouplingMesh *mesh, const double *loc, mcIdType nbOfPoints) const override; + public: + static const char REPR[]; + static constexpr TypeOfField TYPE = ON_NODES_FE; + }; +} diff --git a/src/MEDCoupling/MEDCouplingFieldT.txx b/src/MEDCoupling/MEDCouplingFieldT.txx index 125c314b4..b4c13257f 100644 --- a/src/MEDCoupling/MEDCouplingFieldT.txx +++ b/src/MEDCoupling/MEDCouplingFieldT.txx @@ -449,7 +449,8 @@ namespace MEDCoupling * \ref MEDCoupling::ON_NODES "ON_NODES", * \ref MEDCoupling::ON_GAUSS_PT "ON_GAUSS_PT", * \ref MEDCoupling::ON_GAUSS_NE "ON_GAUSS_NE", - * \ref MEDCoupling::ON_NODES_KR "ON_NODES_KR"). + * \ref MEDCoupling::ON_NODES_KR "ON_NODES_KR", + * \ref MEDCoupling::ON_NODES_FE "ON_NODES_FE"). * * For example, \a this is a field on cells lying on a mesh that have 10 cells, \a partBg contains the following cell ids [3,7,6]. * Then the returned field will lie on mesh having 3 cells and will contain 3 tuples. @@ -512,7 +513,8 @@ namespace MEDCoupling * \ref MEDCoupling::ON_NODES "ON_NODES", * \ref MEDCoupling::ON_GAUSS_PT "ON_GAUSS_PT", * \ref MEDCoupling::ON_GAUSS_NE "ON_GAUSS_NE", - * \ref MEDCoupling::ON_NODES_KR "ON_NODES_KR"). + * \ref MEDCoupling::ON_NODES_KR "ON_NODES_KR", + * \ref MEDCoupling::ON_NODES_FE "ON_NODES_FE"). * * For example, \a this is a field on cells lying on a mesh that have 10 cells, \a part contains following cell ids [3,7,6]. * Then the returned field will lie on mesh having 3 cells and the returned field will contain 3 tuples.
diff --git a/src/MEDCoupling/MEDCouplingNormalizedUnstructuredMesh.hxx b/src/MEDCoupling/MEDCouplingNormalizedUnstructuredMesh.hxx index d3b7ff2a1..7934d0562 100644 --- a/src/MEDCoupling/MEDCouplingNormalizedUnstructuredMesh.hxx +++ b/src/MEDCoupling/MEDCouplingNormalizedUnstructuredMesh.hxx @@ -35,7 +35,7 @@ class MEDCouplingNormalizedUnstructuredMesh public: static const int MY_SPACEDIM=SPACEDIM; static const int MY_MESHDIM=MESHDIM; - typedef mcIdType MyConnType; + using MyConnType = mcIdType; static const INTERP_KERNEL::NumberingPolicy My_numPol=INTERP_KERNEL::ALL_C_MODE; public: MEDCouplingNormalizedUnstructuredMesh(const MEDCoupling::MEDCouplingPointSet *mesh); diff --git a/src/MEDCoupling/MEDCouplingRefCountObject.hxx b/src/MEDCoupling/MEDCouplingRefCountObject.hxx index 674eb9bc5..6dd1ee5c2 100644 --- a/src/MEDCoupling/MEDCouplingRefCountObject.hxx +++ b/src/MEDCoupling/MEDCouplingRefCountObject.hxx @@ -45,7 +45,8 @@ namespace MEDCoupling ON_NODES = 1, ON_GAUSS_PT = 2, ON_GAUSS_NE = 3, - ON_NODES_KR = 4 + ON_NODES_KR = 4, + ON_NODES_FE = 5 } TypeOfField; //! The various temporal discretization of a field diff --git a/src/MEDCoupling_Swig/MEDCouplingCommon.i b/src/MEDCoupling_Swig/MEDCouplingCommon.i index 282847b65..1e25eb35a 100644 --- a/src/MEDCoupling_Swig/MEDCouplingCommon.i +++ b/src/MEDCoupling_Swig/MEDCouplingCommon.i @@ -639,7 +639,8 @@ namespace MEDCoupling ON_NODES = 1, ON_GAUSS_PT = 2, ON_GAUSS_NE = 3, - ON_NODES_KR = 4 + ON_NODES_KR = 4, + ON_NODES_FE = 5 } TypeOfField; typedef enum -- 2.39.2