1 // Copyright (C) 2007-2012 CEA/DEN, EDF R&D
3 // This library is free software; you can redistribute it and/or
4 // modify it under the terms of the GNU Lesser General Public
5 // License as published by the Free Software Foundation; either
6 // version 2.1 of the License.
8 // This library is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
11 // Lesser General Public License for more details.
13 // You should have received a copy of the GNU Lesser General Public
14 // License along with this library; if not, write to the Free Software
15 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
19 // Author : Anthony Geay (CEA/DEN)
21 #include "MEDCouplingRemapper.hxx"
22 #include "MEDCouplingMemArray.hxx"
23 #include "MEDCouplingFieldDouble.hxx"
24 #include "MEDCouplingFieldTemplate.hxx"
25 #include "MEDCouplingFieldDiscretization.hxx"
26 #include "MEDCouplingExtrudedMesh.hxx"
27 #include "MEDCouplingCMesh.hxx"
28 #include "MEDCouplingNormalizedUnstructuredMesh.txx"
29 #include "MEDCouplingNormalizedCartesianMesh.txx"
31 #include "Interpolation1D.txx"
32 #include "Interpolation2DCurve.hxx"
33 #include "Interpolation2D.txx"
34 #include "Interpolation3D.txx"
35 #include "Interpolation3DSurf.hxx"
36 #include "Interpolation2D1D.txx"
37 #include "Interpolation3D2D.txx"
38 #include "InterpolationCU.txx"
39 #include "InterpolationCC.txx"
41 using namespace ParaMEDMEM;
43 MEDCouplingRemapper::MEDCouplingRemapper():_src_mesh(0),_target_mesh(0),_nature_of_deno(NoNature),_time_deno_update(0)
47 MEDCouplingRemapper::~MEDCouplingRemapper()
52 int MEDCouplingRemapper::prepare(const MEDCouplingMesh *srcMesh, const MEDCouplingMesh *targetMesh, const char *method) throw(INTERP_KERNEL::Exception)
55 _src_mesh=const_cast<MEDCouplingMesh *>(srcMesh); _target_mesh=const_cast<MEDCouplingMesh *>(targetMesh);
56 _src_mesh->incrRef(); _target_mesh->incrRef();
57 int meshInterpType=((int)_src_mesh->getType()*16)+(int)_target_mesh->getType();
58 switch(meshInterpType)
60 case 85://Unstructured-Unstructured
61 return prepareUU(method);
62 case 87://Unstructured-Cartesian
63 return prepareUC(method);
64 case 117://Cartesian-Unstructured
65 return prepareCU(method);
66 case 119://Cartesian-Cartesian
67 return prepareCC(method);
68 case 136://Extruded-Extruded
69 return prepareEE(method);
71 throw INTERP_KERNEL::Exception("Not managed type of meshes !");
75 int MEDCouplingRemapper::prepareEx(const MEDCouplingFieldTemplate *src, const MEDCouplingFieldTemplate *target) throw(INTERP_KERNEL::Exception)
77 std::string meth(src->getDiscretization()->getStringRepr());
78 meth+=target->getDiscretization()->getStringRepr();
79 return prepare(src->getMesh(),target->getMesh(),meth.c_str());
83 * This method performs the operation source to target using matrix computed in ParaMEDMEM::MEDCouplingRemapper::prepare method.
84 * If meshes of \b srcField and \b targetField do not match exactly those given into \ref ParaMEDMEM::MEDCouplingRemapper::prepare "prepare method" an exception will be thrown.
86 * \param [in] srcField is the source field from which the interpolation will be done. The mesh into \b srcField should be the same than those specified on ParaMEDMEM::MEDCouplingRemapper::prepare.
87 * \param [out] targetField the destination field with the allocated array in which all tuples will be overwritten.
89 void MEDCouplingRemapper::transfer(const MEDCouplingFieldDouble *srcField, MEDCouplingFieldDouble *targetField, double dftValue) throw(INTERP_KERNEL::Exception)
91 transferUnderground(srcField,targetField,true,dftValue);
95 * This method is equivalent to ParaMEDMEM::MEDCouplingRemapper::transfer except that here \b targetField is a in/out parameter.
96 * If an entity (cell for example) in targetField is not fetched by any entity (cell for example) of \b srcField, the value in targetField is
98 * This method requires that \b targetField was fully defined and allocated. If the array is not allocated an exception will be thrown.
100 * \param [in] srcField is the source field from which the interpolation will be done. The mesh into \b srcField should be the same than those specified on ParaMEDMEM::MEDCouplingRemapper::prepare.
101 * \param [in,out] targetField the destination field with the allocated array in which only tuples whose entities are fetched by interpolation will be overwritten only.
103 void MEDCouplingRemapper::partialTransfer(const MEDCouplingFieldDouble *srcField, MEDCouplingFieldDouble *targetField) throw(INTERP_KERNEL::Exception)
105 transferUnderground(srcField,targetField,false,std::numeric_limits<double>::max());
108 void MEDCouplingRemapper::reverseTransfer(MEDCouplingFieldDouble *srcField, const MEDCouplingFieldDouble *targetField, double dftValue) throw(INTERP_KERNEL::Exception)
110 if(_src_method!=srcField->getDiscretization()->getStringRepr())
111 throw INTERP_KERNEL::Exception("Incoherency with prepare call for source field");
112 if(_target_method!=targetField->getDiscretization()->getStringRepr())
113 throw INTERP_KERNEL::Exception("Incoherency with prepare call for target field");
114 if(srcField->getNature()!=targetField->getNature())
115 throw INTERP_KERNEL::Exception("Natures of fields mismatch !");
116 DataArrayDouble *array=srcField->getArray();
117 int trgNbOfCompo=targetField->getNumberOfComponents();
120 if(trgNbOfCompo!=srcField->getNumberOfComponents())
121 throw INTERP_KERNEL::Exception("Number of components mismatch !");
125 array=DataArrayDouble::New();
126 array->alloc(srcField->getNumberOfTuples(),trgNbOfCompo);
127 srcField->setArray(array);
130 computeDeno(srcField->getNature(),srcField,targetField);
131 double *resPointer=array->getPointer();
132 const double *inputPointer=targetField->getArray()->getConstPointer();
133 computeReverseProduct(inputPointer,trgNbOfCompo,dftValue,resPointer);
136 MEDCouplingFieldDouble *MEDCouplingRemapper::transferField(const MEDCouplingFieldDouble *srcField, double dftValue) throw(INTERP_KERNEL::Exception)
138 if(_src_method!=srcField->getDiscretization()->getStringRepr())
139 throw INTERP_KERNEL::Exception("Incoherency with prepare call for source field");
140 MEDCouplingFieldDouble *ret=MEDCouplingFieldDouble::New(MEDCouplingFieldDiscretization::getTypeOfFieldFromStringRepr(_target_method.c_str()),srcField->getTimeDiscretization());
141 ret->copyAllTinyAttrFrom(srcField);
142 ret->setNature(srcField->getNature());
143 ret->setMesh(_target_mesh);
144 transfer(srcField,ret,dftValue);
148 MEDCouplingFieldDouble *MEDCouplingRemapper::reverseTransferField(const MEDCouplingFieldDouble *targetField, double dftValue) throw(INTERP_KERNEL::Exception)
150 if(_target_method!=targetField->getDiscretization()->getStringRepr())
151 throw INTERP_KERNEL::Exception("Incoherency with prepare call for target field");
152 MEDCouplingFieldDouble *ret=MEDCouplingFieldDouble::New(MEDCouplingFieldDiscretization::getTypeOfFieldFromStringRepr(_src_method.c_str()),targetField->getTimeDiscretization());
153 ret->copyAllTinyAttrFrom(targetField);
154 ret->setNature(targetField->getNature());
155 ret->setMesh(_src_mesh);
156 reverseTransfer(ret,targetField,dftValue);
160 bool MEDCouplingRemapper::setOptionInt(const std::string& key, int value)
162 return INTERP_KERNEL::InterpolationOptions::setOptionInt(key,value);
165 bool MEDCouplingRemapper::setOptionDouble(const std::string& key, double value)
167 return INTERP_KERNEL::InterpolationOptions::setOptionDouble(key,value);
170 bool MEDCouplingRemapper::setOptionString(const std::string& key, const std::string& value)
172 return INTERP_KERNEL::InterpolationOptions::setOptionString(key,value);
175 int MEDCouplingRemapper::prepareUU(const char *method) throw(INTERP_KERNEL::Exception)
177 MEDCouplingUMesh *src_mesh=(MEDCouplingUMesh *)_src_mesh;
178 MEDCouplingUMesh *target_mesh=(MEDCouplingUMesh *)_target_mesh;
179 INTERP_KERNEL::Interpolation<INTERP_KERNEL::Interpolation3D>::checkAndSplitInterpolationMethod(method,_src_method,_target_method);
180 const int srcMeshDim=src_mesh->getMeshDimension();
183 srcSpaceDim=src_mesh->getSpaceDimension();
184 const int trgMeshDim=target_mesh->getMeshDimension();
187 trgSpaceDim=target_mesh->getSpaceDimension();
188 if(trgSpaceDim!=srcSpaceDim)
189 if(trgSpaceDim!=-1 && srcSpaceDim!=-1)
190 throw INTERP_KERNEL::Exception("Incoherent space dimension detected between target and source.");
192 if(srcMeshDim==1 && trgMeshDim==1 && srcSpaceDim==1)
194 MEDCouplingNormalizedUnstructuredMesh<1,1> source_mesh_wrapper(src_mesh);
195 MEDCouplingNormalizedUnstructuredMesh<1,1> target_mesh_wrapper(target_mesh);
196 INTERP_KERNEL::Interpolation1D interpolation(*this);
197 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
199 else if(srcMeshDim==1 && trgMeshDim==1 && srcSpaceDim==2)
201 MEDCouplingNormalizedUnstructuredMesh<2,1> source_mesh_wrapper(src_mesh);
202 MEDCouplingNormalizedUnstructuredMesh<2,1> target_mesh_wrapper(target_mesh);
203 INTERP_KERNEL::Interpolation2DCurve interpolation(*this);
204 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
206 else if(srcMeshDim==2 && trgMeshDim==2 && srcSpaceDim==2)
208 MEDCouplingNormalizedUnstructuredMesh<2,2> source_mesh_wrapper(src_mesh);
209 MEDCouplingNormalizedUnstructuredMesh<2,2> target_mesh_wrapper(target_mesh);
210 INTERP_KERNEL::Interpolation2D interpolation(*this);
211 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
213 else if(srcMeshDim==3 && trgMeshDim==3 && srcSpaceDim==3)
215 MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh);
216 MEDCouplingNormalizedUnstructuredMesh<3,3> target_mesh_wrapper(target_mesh);
217 INTERP_KERNEL::Interpolation3D interpolation(*this);
218 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
220 else if(srcMeshDim==2 && trgMeshDim==2 && srcSpaceDim==3)
222 MEDCouplingNormalizedUnstructuredMesh<3,2> source_mesh_wrapper(src_mesh);
223 MEDCouplingNormalizedUnstructuredMesh<3,2> target_mesh_wrapper(target_mesh);
224 INTERP_KERNEL::Interpolation3DSurf interpolation(*this);
225 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
227 else if(srcMeshDim==3 && trgMeshDim==1 && srcSpaceDim==3)
229 if(getIntersectionType()!=INTERP_KERNEL::PointLocator)
230 throw INTERP_KERNEL::Exception("Invalid interpolation requested between 3D and 1D ! Select PointLocator as intersection type !");
231 MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh);
232 MEDCouplingNormalizedUnstructuredMesh<3,3> target_mesh_wrapper(target_mesh);
233 INTERP_KERNEL::Interpolation3D interpolation(*this);
234 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
236 else if(srcMeshDim==1 && trgMeshDim==3 && srcSpaceDim==3)
238 if(getIntersectionType()!=INTERP_KERNEL::PointLocator)
239 throw INTERP_KERNEL::Exception("Invalid interpolation requested between 3D and 1D ! Select PointLocator as intersection type !");
240 MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh);
241 MEDCouplingNormalizedUnstructuredMesh<3,3> target_mesh_wrapper(target_mesh);
242 INTERP_KERNEL::Interpolation3D interpolation(*this);
243 std::vector<std::map<int,double> > matrixTmp;
244 nbCols=interpolation.interpolateMeshes(target_mesh_wrapper,source_mesh_wrapper,matrixTmp,method);
245 ReverseMatrix(matrixTmp,nbCols,_matrix);
246 nbCols=matrixTmp.size();
248 else if(srcMeshDim==2 && trgMeshDim==1 && srcSpaceDim==2)
250 if(getIntersectionType()==INTERP_KERNEL::PointLocator)
252 MEDCouplingNormalizedUnstructuredMesh<2,2> source_mesh_wrapper(src_mesh);
253 MEDCouplingNormalizedUnstructuredMesh<2,2> target_mesh_wrapper(target_mesh);
254 INTERP_KERNEL::Interpolation2D interpolation(*this);
255 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
259 MEDCouplingNormalizedUnstructuredMesh<2,2> source_mesh_wrapper(src_mesh);
260 MEDCouplingNormalizedUnstructuredMesh<2,2> target_mesh_wrapper(target_mesh);
261 INTERP_KERNEL::Interpolation2D1D interpolation(*this);
262 std::vector<std::map<int,double> > matrixTmp;
263 nbCols=interpolation.interpolateMeshes(target_mesh_wrapper,source_mesh_wrapper,matrixTmp,method);
264 ReverseMatrix(matrixTmp,nbCols,_matrix);
265 nbCols=matrixTmp.size();
266 INTERP_KERNEL::Interpolation2D1D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces();
267 if(!duplicateFaces.empty())
269 std::ostringstream oss; oss << "An unexpected situation happend ! For the following 1D Cells are part of edges shared by 2D cells :\n";
270 for(std::map<int,std::set<int> >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++)
272 oss << "1D Cell #" << (*it).first << " is part of common edge of following 2D cells ids : ";
273 std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator<int>(oss," "));
279 else if(srcMeshDim==1 && trgMeshDim==2 && srcSpaceDim==2)
281 if(getIntersectionType()==INTERP_KERNEL::PointLocator)
283 MEDCouplingNormalizedUnstructuredMesh<2,2> source_mesh_wrapper(src_mesh);
284 MEDCouplingNormalizedUnstructuredMesh<2,2> target_mesh_wrapper(target_mesh);
285 INTERP_KERNEL::Interpolation2D interpolation(*this);
286 std::vector<std::map<int,double> > matrixTmp;
287 nbCols=interpolation.interpolateMeshes(target_mesh_wrapper,source_mesh_wrapper,matrixTmp,method);
288 ReverseMatrix(matrixTmp,nbCols,_matrix);
289 nbCols=matrixTmp.size();
293 MEDCouplingNormalizedUnstructuredMesh<2,2> source_mesh_wrapper(src_mesh);
294 MEDCouplingNormalizedUnstructuredMesh<2,2> target_mesh_wrapper(target_mesh);
295 INTERP_KERNEL::Interpolation2D1D interpolation(*this);
296 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
297 INTERP_KERNEL::Interpolation2D1D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces();
298 if(!duplicateFaces.empty())
300 std::ostringstream oss; oss << "An unexpected situation happend ! For the following 1D Cells are part of edges shared by 2D cells :\n";
301 for(std::map<int,std::set<int> >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++)
303 oss << "1D Cell #" << (*it).first << " is part of common edge of following 2D cells ids : ";
304 std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator<int>(oss," "));
310 else if(srcMeshDim==2 && trgMeshDim==3 && srcSpaceDim==3)
312 MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh);
313 MEDCouplingNormalizedUnstructuredMesh<3,3> target_mesh_wrapper(target_mesh);
314 INTERP_KERNEL::Interpolation3D2D interpolation(*this);
315 nbCols=interpolation.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,_matrix,method);
316 INTERP_KERNEL::Interpolation3D2D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces();
317 if(!duplicateFaces.empty())
319 std::ostringstream oss; oss << "An unexpected situation happend ! For the following 2D Cells are part of edges shared by 3D cells :\n";
320 for(std::map<int,std::set<int> >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++)
322 oss << "2D Cell #" << (*it).first << " is part of common face of following 3D cells ids : ";
323 std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator<int>(oss," "));
328 else if(srcMeshDim==3 && trgMeshDim==2 && srcSpaceDim==3)
330 MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh);
331 MEDCouplingNormalizedUnstructuredMesh<3,3> target_mesh_wrapper(target_mesh);
332 INTERP_KERNEL::Interpolation3D2D interpolation(*this);
333 std::vector<std::map<int,double> > matrixTmp;
334 nbCols=interpolation.interpolateMeshes(target_mesh_wrapper,source_mesh_wrapper,matrixTmp,method);
335 ReverseMatrix(matrixTmp,nbCols,_matrix);
336 nbCols=matrixTmp.size();
337 INTERP_KERNEL::Interpolation3D2D::DuplicateFacesType duplicateFaces=interpolation.retrieveDuplicateFaces();
338 if(!duplicateFaces.empty())
340 std::ostringstream oss; oss << "An unexpected situation happend ! For the following 2D Cells are part of edges shared by 3D cells :\n";
341 for(std::map<int,std::set<int> >::const_iterator it=duplicateFaces.begin();it!=duplicateFaces.end();it++)
343 oss << "2D Cell #" << (*it).first << " is part of common face of following 3D cells ids : ";
344 std::copy((*it).second.begin(),(*it).second.end(),std::ostream_iterator<int>(oss," "));
349 else if(trgMeshDim==-1)
351 if(srcMeshDim==2 && srcSpaceDim==2)
353 MEDCouplingNormalizedUnstructuredMesh<2,2> source_mesh_wrapper(src_mesh);
354 INTERP_KERNEL::Interpolation2D interpolation(*this);
355 nbCols=interpolation.toIntegralUniform(source_mesh_wrapper,_matrix,_src_method.c_str());
357 else if(srcMeshDim==3 && srcSpaceDim==3)
359 MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(src_mesh);
360 INTERP_KERNEL::Interpolation3D interpolation(*this);
361 nbCols=interpolation.toIntegralUniform(source_mesh_wrapper,_matrix,_src_method.c_str());
363 else if(srcMeshDim==2 && srcSpaceDim==3)
365 MEDCouplingNormalizedUnstructuredMesh<3,2> source_mesh_wrapper(src_mesh);
366 INTERP_KERNEL::Interpolation3DSurf interpolation(*this);
367 nbCols=interpolation.toIntegralUniform(source_mesh_wrapper,_matrix,_src_method.c_str());
370 throw INTERP_KERNEL::Exception("No interpolation available for the given mesh and space dimension of source mesh to -1D targetMesh");
372 else if(srcMeshDim==-1)
374 if(trgMeshDim==2 && trgSpaceDim==2)
376 MEDCouplingNormalizedUnstructuredMesh<2,2> source_mesh_wrapper(target_mesh);
377 INTERP_KERNEL::Interpolation2D interpolation(*this);
378 nbCols=interpolation.fromIntegralUniform(source_mesh_wrapper,_matrix,_target_method.c_str());
380 else if(trgMeshDim==3 && trgSpaceDim==3)
382 MEDCouplingNormalizedUnstructuredMesh<3,3> source_mesh_wrapper(target_mesh);
383 INTERP_KERNEL::Interpolation3D interpolation(*this);
384 nbCols=interpolation.fromIntegralUniform(source_mesh_wrapper,_matrix,_target_method.c_str());
386 else if(trgMeshDim==2 && trgSpaceDim==3)
388 MEDCouplingNormalizedUnstructuredMesh<3,2> source_mesh_wrapper(target_mesh);
389 INTERP_KERNEL::Interpolation3DSurf interpolation(*this);
390 nbCols=interpolation.fromIntegralUniform(source_mesh_wrapper,_matrix,_target_method.c_str());
393 throw INTERP_KERNEL::Exception("No interpolation available for the given mesh and space dimension of source mesh from -1D sourceMesh");
396 throw INTERP_KERNEL::Exception("No interpolation available for the given mesh and space dimension");
397 _deno_multiply.clear();
398 _deno_multiply.resize(_matrix.size());
399 _deno_reverse_multiply.clear();
400 _deno_reverse_multiply.resize(nbCols);
405 int MEDCouplingRemapper::prepareEE(const char *method) throw(INTERP_KERNEL::Exception)
407 MEDCouplingExtrudedMesh *src_mesh=(MEDCouplingExtrudedMesh *)_src_mesh;
408 MEDCouplingExtrudedMesh *target_mesh=(MEDCouplingExtrudedMesh *)_target_mesh;
409 std::string methC(method);
411 throw INTERP_KERNEL::Exception("Only P0P0 method implemented for Extruded/Extruded meshes !");
412 INTERP_KERNEL::Interpolation<INTERP_KERNEL::Interpolation3D>::checkAndSplitInterpolationMethod(method,_src_method,_target_method);
413 MEDCouplingNormalizedUnstructuredMesh<3,2> source_mesh_wrapper(src_mesh->getMesh2D());
414 MEDCouplingNormalizedUnstructuredMesh<3,2> target_mesh_wrapper(target_mesh->getMesh2D());
415 INTERP_KERNEL::Interpolation3DSurf interpolation2D(*this);
416 std::vector<std::map<int,double> > matrix2D;
417 int nbCols2D=interpolation2D.interpolateMeshes(source_mesh_wrapper,target_mesh_wrapper,matrix2D,method);
418 MEDCouplingUMesh *s1D,*t1D;
420 MEDCouplingExtrudedMesh::Project1DMeshes(src_mesh->getMesh1D(),target_mesh->getMesh1D(),getPrecision(),s1D,t1D,v);
421 MEDCouplingNormalizedUnstructuredMesh<1,1> s1DWrapper(s1D);
422 MEDCouplingNormalizedUnstructuredMesh<1,1> t1DWrapper(t1D);
423 std::vector<std::map<int,double> > matrix1D;
424 INTERP_KERNEL::Interpolation1D interpolation1D(*this);
425 int nbCols1D=interpolation1D.interpolateMeshes(s1DWrapper,t1DWrapper,matrix1D,method);
428 buildFinalInterpolationMatrixByConvolution(matrix1D,matrix2D,src_mesh->getMesh3DIds()->getConstPointer(),nbCols2D,nbCols1D,
429 target_mesh->getMesh3DIds()->getConstPointer());
431 _deno_multiply.clear();
432 _deno_multiply.resize(_matrix.size());
433 _deno_reverse_multiply.clear();
434 _deno_reverse_multiply.resize(nbCols2D*nbCols1D);
439 int MEDCouplingRemapper::prepareUC(const char *method) throw(INTERP_KERNEL::Exception)
441 std::string methodCpp(method);
442 if(methodCpp!="P0P0")
443 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareUC : only P0P0 interpolation supported for the moment !");
444 INTERP_KERNEL::Interpolation<INTERP_KERNEL::Interpolation3D>::checkAndSplitInterpolationMethod(method,_src_method,_target_method);
445 MEDCouplingUMesh *src_mesh=static_cast<MEDCouplingUMesh *>(_src_mesh);
446 MEDCouplingCMesh *target_mesh=static_cast<MEDCouplingCMesh *>(_target_mesh);
447 const int srcMeshDim=src_mesh->getMeshDimension();
448 const int srcSpceDim=src_mesh->getSpaceDimension();
449 const int trgMeshDim=target_mesh->getMeshDimension();
450 if(srcMeshDim!=srcSpceDim || srcMeshDim!=trgMeshDim)
451 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareUC : space dim of src unstructured should be equal to mesh dim of src unstructured and should be equal also equal to trg cartesian dimension !");
452 std::vector<std::map<int,double> > res;
457 MEDCouplingNormalizedCartesianMesh<1> targetWrapper(target_mesh);
458 MEDCouplingNormalizedUnstructuredMesh<1,1> sourceWrapper(src_mesh);
459 INTERP_KERNEL::InterpolationCU myInterpolator(*this);
460 myInterpolator.interpolateMeshes(targetWrapper,sourceWrapper,res,"P0P0");
465 MEDCouplingNormalizedCartesianMesh<2> targetWrapper(target_mesh);
466 MEDCouplingNormalizedUnstructuredMesh<2,2> sourceWrapper(src_mesh);
467 INTERP_KERNEL::InterpolationCU myInterpolator(*this);
468 myInterpolator.interpolateMeshes(targetWrapper,sourceWrapper,res,"P0P0");
473 MEDCouplingNormalizedCartesianMesh<3> targetWrapper(target_mesh);
474 MEDCouplingNormalizedUnstructuredMesh<3,3> sourceWrapper(src_mesh);
475 INTERP_KERNEL::InterpolationCU myInterpolator(*this);
476 myInterpolator.interpolateMeshes(targetWrapper,sourceWrapper,res,"P0P0");
480 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareUC : only dimension 1 2 or 3 supported !");
482 ReverseMatrix(res,target_mesh->getNumberOfCells(),_matrix);
483 nullifiedTinyCoeffInCrudeMatrixAbs(0.);
485 _deno_multiply.clear();
486 _deno_multiply.resize(_matrix.size());
487 _deno_reverse_multiply.clear();
488 _deno_reverse_multiply.resize(src_mesh->getNumberOfCells());
493 int MEDCouplingRemapper::prepareCU(const char *method) throw(INTERP_KERNEL::Exception)
495 std::string methodCpp(method);
496 if(methodCpp!="P0P0")
497 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareCU : only P0P0 interpolation supported for the moment !");
498 INTERP_KERNEL::Interpolation<INTERP_KERNEL::Interpolation3D>::checkAndSplitInterpolationMethod(method,_src_method,_target_method);
499 MEDCouplingCMesh *src_mesh=static_cast<MEDCouplingCMesh *>(_src_mesh);
500 MEDCouplingUMesh *target_mesh=static_cast<MEDCouplingUMesh *>(_target_mesh);
501 const int srcMeshDim=src_mesh->getMeshDimension();
502 const int trgMeshDim=target_mesh->getMeshDimension();
503 const int trgSpceDim=target_mesh->getSpaceDimension();
504 if(trgMeshDim!=trgSpceDim || trgMeshDim!=srcMeshDim)
505 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareCU : space dim of target unstructured should be equal to mesh dim of target unstructured and should be equal also equal to source cartesian dimension !");
510 MEDCouplingNormalizedCartesianMesh<1> sourceWrapper(src_mesh);
511 MEDCouplingNormalizedUnstructuredMesh<1,1> targetWrapper(target_mesh);
512 INTERP_KERNEL::InterpolationCU myInterpolator(*this);
513 myInterpolator.interpolateMeshes(sourceWrapper,targetWrapper,_matrix,"P0P0");
518 MEDCouplingNormalizedCartesianMesh<2> sourceWrapper(src_mesh);
519 MEDCouplingNormalizedUnstructuredMesh<2,2> targetWrapper(target_mesh);
520 INTERP_KERNEL::InterpolationCU myInterpolator(*this);
521 myInterpolator.interpolateMeshes(sourceWrapper,targetWrapper,_matrix,"P0P0");
526 MEDCouplingNormalizedCartesianMesh<3> sourceWrapper(src_mesh);
527 MEDCouplingNormalizedUnstructuredMesh<3,3> targetWrapper(target_mesh);
528 INTERP_KERNEL::InterpolationCU myInterpolator(*this);
529 myInterpolator.interpolateMeshes(sourceWrapper,targetWrapper,_matrix,"P0P0");
533 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareCU : only dimension 1 2 or 3 supported !");
535 nullifiedTinyCoeffInCrudeMatrixAbs(0.);
537 _deno_multiply.clear();
538 _deno_multiply.resize(_matrix.size());
539 _deno_reverse_multiply.clear();
540 _deno_reverse_multiply.resize(src_mesh->getNumberOfCells());
545 int MEDCouplingRemapper::prepareCC(const char *method) throw(INTERP_KERNEL::Exception)
547 std::string methodCpp(method);
548 if(methodCpp!="P0P0")
549 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareCC : only P0P0 interpolation supported for the moment !");
550 INTERP_KERNEL::Interpolation<INTERP_KERNEL::Interpolation3D>::checkAndSplitInterpolationMethod(method,_src_method,_target_method);
551 MEDCouplingCMesh *src_mesh=static_cast<MEDCouplingCMesh *>(_src_mesh);
552 MEDCouplingCMesh *target_mesh=static_cast<MEDCouplingCMesh *>(_target_mesh);
553 const int srcMeshDim=src_mesh->getMeshDimension();
554 const int trgMeshDim=target_mesh->getMeshDimension();
555 if(trgMeshDim!=srcMeshDim)
556 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareCC : dim of target cartesian should be equal to dim of source cartesian dimension !");
561 MEDCouplingNormalizedCartesianMesh<1> sourceWrapper(src_mesh);
562 MEDCouplingNormalizedCartesianMesh<1> targetWrapper(target_mesh);
563 INTERP_KERNEL::InterpolationCC myInterpolator(*this);
564 myInterpolator.interpolateMeshes(sourceWrapper,targetWrapper,_matrix,"P0P0");
569 MEDCouplingNormalizedCartesianMesh<2> sourceWrapper(src_mesh);
570 MEDCouplingNormalizedCartesianMesh<2> targetWrapper(target_mesh);
571 INTERP_KERNEL::InterpolationCC myInterpolator(*this);
572 myInterpolator.interpolateMeshes(sourceWrapper,targetWrapper,_matrix,"P0P0");
577 MEDCouplingNormalizedCartesianMesh<3> sourceWrapper(src_mesh);
578 MEDCouplingNormalizedCartesianMesh<3> targetWrapper(target_mesh);
579 INTERP_KERNEL::InterpolationCC myInterpolator(*this);
580 myInterpolator.interpolateMeshes(sourceWrapper,targetWrapper,_matrix,"P0P0");
584 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::prepareCC : only dimension 1 2 or 3 supported !");
586 nullifiedTinyCoeffInCrudeMatrixAbs(0.);
588 _deno_multiply.clear();
589 _deno_multiply.resize(_matrix.size());
590 _deno_reverse_multiply.clear();
591 _deno_reverse_multiply.resize(src_mesh->getNumberOfCells());
596 void MEDCouplingRemapper::updateTime() const
600 void MEDCouplingRemapper::releaseData(bool matrixSuppression)
603 _src_mesh->decrRef();
605 _target_mesh->decrRef();
608 if(matrixSuppression)
611 _deno_multiply.clear();
612 _deno_reverse_multiply.clear();
616 void MEDCouplingRemapper::transferUnderground(const MEDCouplingFieldDouble *srcField, MEDCouplingFieldDouble *targetField, bool isDftVal, double dftValue) throw(INTERP_KERNEL::Exception)
618 if(_src_method!=srcField->getDiscretization()->getStringRepr())
619 throw INTERP_KERNEL::Exception("Incoherency with prepare call for source field");
620 if(_target_method!=targetField->getDiscretization()->getStringRepr())
621 throw INTERP_KERNEL::Exception("Incoherency with prepare call for target field");
622 if(srcField->getNature()!=targetField->getNature())
623 throw INTERP_KERNEL::Exception("Natures of fields mismatch !");
624 DataArrayDouble *array=targetField->getArray();
625 int srcNbOfCompo=srcField->getNumberOfComponents();
628 if(srcNbOfCompo!=targetField->getNumberOfComponents())
629 throw INTERP_KERNEL::Exception("Number of components mismatch !");
634 throw INTERP_KERNEL::Exception("MEDCouplingRemapper::partialTransfer : This method requires that the array of target field exists ! Allocate it or call MEDCouplingRemapper::transfer instead !");
635 array=DataArrayDouble::New();
636 array->alloc(targetField->getNumberOfTuples(),srcNbOfCompo);
637 targetField->setArray(array);
640 computeDeno(srcField->getNature(),srcField,targetField);
641 double *resPointer=array->getPointer();
642 const double *inputPointer=srcField->getArray()->getConstPointer();
643 computeProduct(inputPointer,srcNbOfCompo,isDftVal,dftValue,resPointer);
646 void MEDCouplingRemapper::computeDeno(NatureOfField nat, const MEDCouplingFieldDouble *srcField, const MEDCouplingFieldDouble *trgField)
649 return computeDenoFromScratch(nat,srcField,trgField);
650 else if(nat!=_nature_of_deno)
651 return computeDenoFromScratch(nat,srcField,trgField);
652 else if(nat==_nature_of_deno && _time_deno_update!=getTimeOfThis())
653 return computeDenoFromScratch(nat,srcField,trgField);
656 void MEDCouplingRemapper::computeDenoFromScratch(NatureOfField nat, const MEDCouplingFieldDouble *srcField, const MEDCouplingFieldDouble *trgField) throw(INTERP_KERNEL::Exception)
659 _time_deno_update=getTimeOfThis();
660 switch(_nature_of_deno)
662 case ConservativeVolumic:
664 ComputeRowSumAndColSum(_matrix,_deno_multiply,_deno_reverse_multiply);
669 MEDCouplingFieldDouble *deno=srcField->getDiscretization()->getMeasureField(srcField->getMesh(),getMeasureAbsStatus());
670 MEDCouplingFieldDouble *denoR=trgField->getDiscretization()->getMeasureField(trgField->getMesh(),getMeasureAbsStatus());
671 const double *denoPtr=deno->getArray()->getConstPointer();
672 const double *denoRPtr=denoR->getArray()->getConstPointer();
673 if(trgField->getMesh()->getMeshDimension()==-1)
675 double *denoRPtr2=denoR->getArray()->getPointer();
676 denoRPtr2[0]=std::accumulate(denoPtr,denoPtr+deno->getNumberOfTuples(),0.);
678 if(srcField->getMesh()->getMeshDimension()==-1)
680 double *denoPtr2=deno->getArray()->getPointer();
681 denoPtr2[0]=std::accumulate(denoRPtr,denoRPtr+denoR->getNumberOfTuples(),0.);
684 for(std::vector<std::map<int,double> >::const_iterator iter1=_matrix.begin();iter1!=_matrix.end();iter1++,idx++)
685 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
687 _deno_multiply[idx][(*iter2).first]=denoPtr[(*iter2).first];
688 _deno_reverse_multiply[(*iter2).first][idx]=denoRPtr[idx];
694 case IntegralGlobConstraint:
696 ComputeColSumAndRowSum(_matrix,_deno_multiply,_deno_reverse_multiply);
701 MEDCouplingFieldDouble *deno=trgField->getDiscretization()->getMeasureField(trgField->getMesh(),getMeasureAbsStatus());
702 MEDCouplingFieldDouble *denoR=srcField->getDiscretization()->getMeasureField(srcField->getMesh(),getMeasureAbsStatus());
703 const double *denoPtr=deno->getArray()->getConstPointer();
704 const double *denoRPtr=denoR->getArray()->getConstPointer();
705 if(trgField->getMesh()->getMeshDimension()==-1)
707 double *denoRPtr2=denoR->getArray()->getPointer();
708 denoRPtr2[0]=std::accumulate(denoPtr,denoPtr+deno->getNumberOfTuples(),0.);
710 if(srcField->getMesh()->getMeshDimension()==-1)
712 double *denoPtr2=deno->getArray()->getPointer();
713 denoPtr2[0]=std::accumulate(denoRPtr,denoRPtr+denoR->getNumberOfTuples(),0.);
716 for(std::vector<std::map<int,double> >::const_iterator iter1=_matrix.begin();iter1!=_matrix.end();iter1++,idx++)
717 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
719 _deno_multiply[idx][(*iter2).first]=denoPtr[idx];
720 _deno_reverse_multiply[(*iter2).first][idx]=denoRPtr[(*iter2).first];
727 throw INTERP_KERNEL::Exception("No nature specified ! Select one !");
731 void MEDCouplingRemapper::computeProduct(const double *inputPointer, int inputNbOfCompo, bool isDftVal, double dftValue, double *resPointer)
734 double *tmp=new double[inputNbOfCompo];
735 for(std::vector<std::map<int,double> >::const_iterator iter1=_matrix.begin();iter1!=_matrix.end();iter1++,idx++)
740 std::fill(resPointer+idx*inputNbOfCompo,resPointer+(idx+1)*inputNbOfCompo,dftValue);
744 std::fill(resPointer+idx*inputNbOfCompo,resPointer+(idx+1)*inputNbOfCompo,0.);
745 std::map<int,double>::const_iterator iter3=_deno_multiply[idx].begin();
746 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++,iter3++)
748 std::transform(inputPointer+(*iter2).first*inputNbOfCompo,inputPointer+((*iter2).first+1)*inputNbOfCompo,tmp,std::bind2nd(std::multiplies<double>(),(*iter2).second/(*iter3).second));
749 std::transform(tmp,tmp+inputNbOfCompo,resPointer+idx*inputNbOfCompo,resPointer+idx*inputNbOfCompo,std::plus<double>());
755 void MEDCouplingRemapper::computeReverseProduct(const double *inputPointer, int inputNbOfCompo, double dftValue, double *resPointer)
757 std::vector<bool> isReached(_deno_reverse_multiply.size(),false);
759 double *tmp=new double[inputNbOfCompo];
760 std::fill(resPointer,resPointer+inputNbOfCompo*_deno_reverse_multiply.size(),0.);
761 for(std::vector<std::map<int,double> >::const_iterator iter1=_matrix.begin();iter1!=_matrix.end();iter1++,idx++)
763 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
765 isReached[(*iter2).first]=true;
766 std::transform(inputPointer+idx*inputNbOfCompo,inputPointer+(idx+1)*inputNbOfCompo,tmp,std::bind2nd(std::multiplies<double>(),(*iter2).second/_deno_reverse_multiply[(*iter2).first][idx]));
767 std::transform(tmp,tmp+inputNbOfCompo,resPointer+((*iter2).first)*inputNbOfCompo,resPointer+((*iter2).first)*inputNbOfCompo,std::plus<double>());
772 for(std::vector<bool>::const_iterator iter3=isReached.begin();iter3!=isReached.end();iter3++,idx++)
774 std::fill(resPointer+idx*inputNbOfCompo,resPointer+(idx+1)*inputNbOfCompo,dftValue);
777 void MEDCouplingRemapper::ReverseMatrix(const std::vector<std::map<int,double> >& matIn, int nbColsMatIn, std::vector<std::map<int,double> >& matOut)
779 matOut.resize(nbColsMatIn);
781 for(std::vector<std::map<int,double> >::const_iterator iter1=matIn.begin();iter1!=matIn.end();iter1++,id++)
782 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
783 matOut[(*iter2).first][id]=(*iter2).second;
786 void MEDCouplingRemapper::ComputeRowSumAndColSum(const std::vector<std::map<int,double> >& matrixDeno,
787 std::vector<std::map<int,double> >& deno, std::vector<std::map<int,double> >& denoReverse)
789 std::map<int,double> values;
791 for(std::vector<std::map<int,double> >::const_iterator iter1=matrixDeno.begin();iter1!=matrixDeno.end();iter1++,idx++)
794 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
796 sum+=(*iter2).second;
797 values[(*iter2).first]+=(*iter2).second;
799 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
800 deno[idx][(*iter2).first]=sum;
803 for(std::vector<std::map<int,double> >::const_iterator iter1=matrixDeno.begin();iter1!=matrixDeno.end();iter1++,idx++)
805 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
806 denoReverse[(*iter2).first][idx]=values[(*iter2).first];
810 void MEDCouplingRemapper::ComputeColSumAndRowSum(const std::vector<std::map<int,double> >& matrixDeno,
811 std::vector<std::map<int,double> >& deno, std::vector<std::map<int,double> >& denoReverse)
813 std::map<int,double> values;
815 for(std::vector<std::map<int,double> >::const_iterator iter1=matrixDeno.begin();iter1!=matrixDeno.end();iter1++,idx++)
818 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
820 sum+=(*iter2).second;
821 values[(*iter2).first]+=(*iter2).second;
823 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
824 denoReverse[(*iter2).first][idx]=sum;
827 for(std::vector<std::map<int,double> >::const_iterator iter1=matrixDeno.begin();iter1!=matrixDeno.end();iter1++,idx++)
829 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
830 deno[idx][(*iter2).first]=values[(*iter2).first];
834 void MEDCouplingRemapper::buildFinalInterpolationMatrixByConvolution(const std::vector< std::map<int,double> >& m1D,
835 const std::vector< std::map<int,double> >& m2D,
836 const int *corrCellIdSrc, int nbOf2DCellsSrc, int nbOf1DCellsSrc,
837 const int *corrCellIdTrg)
839 int nbOf2DCellsTrg=m2D.size();
840 int nbOf1DCellsTrg=m1D.size();
841 int nbOf3DCellsTrg=nbOf2DCellsTrg*nbOf1DCellsTrg;
842 _matrix.resize(nbOf3DCellsTrg);
844 for(std::vector< std::map<int,double> >::const_iterator iter2R=m2D.begin();iter2R!=m2D.end();iter2R++,id2R++)
846 for(std::map<int,double>::const_iterator iter2C=(*iter2R).begin();iter2C!=(*iter2R).end();iter2C++)
849 for(std::vector< std::map<int,double> >::const_iterator iter1R=m1D.begin();iter1R!=m1D.end();iter1R++,id1R++)
851 for(std::map<int,double>::const_iterator iter1C=(*iter1R).begin();iter1C!=(*iter1R).end();iter1C++)
853 _matrix[corrCellIdTrg[id1R*nbOf2DCellsTrg+id2R]][corrCellIdSrc[(*iter1C).first*nbOf2DCellsSrc+(*iter2C).first]]=(*iter1C).second*((*iter2C).second);
860 void MEDCouplingRemapper::PrintMatrix(const std::vector<std::map<int,double> >& m)
863 for(std::vector<std::map<int,double> >::const_iterator iter1=m.begin();iter1!=m.end();iter1++,id++)
865 std::cout << "Target Cell # " << id << " : ";
866 for(std::map<int,double>::const_iterator iter2=(*iter1).begin();iter2!=(*iter1).end();iter2++)
867 std::cout << "(" << (*iter2).first << "," << (*iter2).second << "), ";
868 std::cout << std::endl;
872 const std::vector<std::map<int,double> >& MEDCouplingRemapper::getCrudeMatrix() const
878 * This method is supposed to be called , if needed, right after MEDCouplingRemapper::prepare or MEDCouplingRemapper::prepareEx.
879 * If not the behaviour is unpredictable.
880 * This method works on precomputed \a this->_matrix. All coefficients in the matrix is lower than \a maxValAbs this coefficient is
881 * set to 0. That is to say that its entry disappear from the map storing the corresponding row in the data storage of sparse crude matrix.
882 * This method is useful to correct at a high level some problems linked to precision. Indeed, with some \ref NatureOfField "natures of field" some threshold effect
885 * \param [in] maxValAbs is a limit behind which a coefficient is set to 0. \a maxValAbs is expected to be positive, if not this method do nothing.
886 * \return a positive value that tells the number of coefficients put to 0. The 0 returned value means that the matrix has remained unchanged.
887 * \sa MEDCouplingRemapper::nullifiedTinyCoeffInCrudeMatrix
889 int MEDCouplingRemapper::nullifiedTinyCoeffInCrudeMatrixAbs(double maxValAbs) throw(INTERP_KERNEL::Exception)
892 std::vector<std::map<int,double> > matrixNew(_matrix.size());
894 for(std::vector<std::map<int,double> >::const_iterator it1=_matrix.begin();it1!=_matrix.end();it1++,i++)
896 std::map<int,double>& rowNew=matrixNew[i];
897 for(std::map<int,double>::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
899 if(fabs((*it2).second)>maxValAbs)
900 rowNew[(*it2).first]=(*it2).second;
911 * This method is supposed to be called , if needed, right after MEDCouplingRemapper::prepare or MEDCouplingRemapper::prepareEx.
912 * If not the behaviour is unpredictable.
913 * This method works on precomputed \a this->_matrix. All coefficients in the matrix is lower than delta multiplied by \a scaleFactor this coefficient is
914 * set to 0. That is to say that its entry disappear from the map storing the corresponding row in the data storage of sparse crude matrix.
915 * delta is the value returned by MEDCouplingRemapper::getMaxValueInCrudeMatrix method.
916 * This method is useful to correct at a high level some problems linked to precision. Indeed, with some \ref NatureOfField "natures of field" some threshold effect
919 * \param [in] scaleFactor is the scale factor from which coefficients lower than \a scaleFactor times range width of coefficients are set to zero.
920 * \return a positive value that tells the number of coefficients put to 0. The 0 returned value means that the matrix has remained unchanged. If -1 is returned it means
921 * that all coefficients are null.
922 * \sa MEDCouplingRemapper::nullifiedTinyCoeffInCrudeMatrixAbs
924 int MEDCouplingRemapper::nullifiedTinyCoeffInCrudeMatrix(double scaleFactor) throw(INTERP_KERNEL::Exception)
926 double maxVal=getMaxValueInCrudeMatrix();
929 return nullifiedTinyCoeffInCrudeMatrixAbs(scaleFactor*maxVal);
933 * This method is supposed to be called , if needed, right after MEDCouplingRemapper::prepare or MEDCouplingRemapper::prepareEx.
934 * If not the behaviour is unpredictable.
935 * This method returns the maximum of the absolute values of coefficients into the sparse crude matrix.
936 * The returned value is positive.
938 double MEDCouplingRemapper::getMaxValueInCrudeMatrix() const throw(INTERP_KERNEL::Exception)
941 for(std::vector<std::map<int,double> >::const_iterator it1=_matrix.begin();it1!=_matrix.end();it1++)
942 for(std::map<int,double>::const_iterator it2=(*it1).begin();it2!=(*it1).end();it2++)
943 if(fabs((*it2).second)>ret)
944 ret=fabs((*it2).second);