1 // Copyright (C) 2021 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, or (at your option) any later version.
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 (EDF R&D)
21 #ifndef vtkSedimentDeposit_h__
22 #define vtkSedimentDeposit_h__
24 #include <vtkDataObjectAlgorithm.h>
25 #include <vtkUnstructuredGrid.h>
26 #include <vtkSetGet.h>
27 #include <vtkDataSet.h>
29 #include <vtkAdjacentVertexIterator.h>
30 #include <vtkAlgorithmOutput.h>
32 #include <vtkCellData.h>
33 #include <vtkCellType.h>
34 #include <vtkCharArray.h>
35 #include <vtkDataArraySelection.h>
36 #include <vtkDataObjectTreeIterator.h>
37 #include <vtkDataSet.h>
38 #include <vtkDataSetAttributes.h>
39 #include <vtkDemandDrivenPipeline.h>
40 #include <vtkDoubleArray.h>
41 #include <vtkExecutive.h>
42 #include <vtkDoubleArray.h>
43 #include <vtkInEdgeIterator.h>
44 #include <vtkInformation.h>
45 #include <vtkInformationDataObjectKey.h>
46 #include <vtkInformationStringKey.h>
47 #include <vtkInformationVector.h>
48 #include <vtkIntArray.h>
49 #include <vtkMultiBlockDataSet.h>
50 #include <vtkMutableDirectedGraph.h>
52 #include <vtkObjectFactory.h>
53 #include <vtkPointData.h>
54 #include <vtkPolyData.h>
55 #include <vtkResampleWithDataSet.h>
56 #include <vtkStreamingDemandDrivenPipeline.h>
57 #include <vtkStringArray.h>
59 #include <vtkTimeStamp.h>
60 #include <vtkUnsignedCharArray.h>
61 #include <vtkUnstructuredGrid.h>
62 #include <vtkVariantArray.h>
63 #include <vtkWarpScalar.h>
65 #include "VTKToMEDMem.h"
73 class VTK_EXPORT vtkSedimentDeposit : public vtkDataObjectAlgorithm
76 static vtkSedimentDeposit *New();
77 vtkTypeMacro(vtkSedimentDeposit, vtkDataObjectAlgorithm);
78 void PrintSelf(ostream &os, vtkIndent indent) override;
80 void SetSourceData(vtkDataObject *input);
81 void SetSourceConnection(vtkAlgorithmOutput *algOutput);
83 int FillOutputPortInformation(int, vtkInformation *) override;
87 ~vtkSedimentDeposit() override = default;
89 int RequestInformation(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override;
90 int RequestData(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override;
91 int RequestUpdateExtent(vtkInformation *, vtkInformationVector **, vtkInformationVector *) override;
93 int NumberOfTimeSteps;
96 class VTK_EXPORT vtkInternal
99 vtkInternal(int curveId, int nbOfCurves) :_curve_id(curveId), _nb_of_curves(nbOfCurves) {}
100 void pushData(double timeStep, double positiveValue, double negativeValue) { _data.emplace_back(timeStep, positiveValue, negativeValue); }
101 void fillTable(vtkTable *table) const;
102 void analyzeInputDataSets(vtkUnstructuredGrid *ds1, vtkDataSet *ds2);
103 bool computationNeeded() const;
104 MEDCoupling::MCAuto<MEDCoupling::MEDCouplingUMesh> &meshOrigin() { return _meshorig; }
105 MEDCoupling::MCAuto<MEDCoupling::DataArrayIdType> &untouched2DCells() { return _untouched_2d_cells; }
106 MEDCoupling::MCAuto<MEDCoupling::DataArrayIdType> &cellsAtBoundary() { return _cells_at_boundary_origin; }
107 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> ¢ers() { return _centers; }
108 MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> &measure() { return _vol; }
109 std::string getReprDependingPos(const std::string& origName) const;
112 std::vector<std::tuple<double, double, double>> _data;
113 vtkMTimeType _mt1 = 0;
114 vtkMTimeType _mt2 = 0;
117 bool _recomputationOfMatrixNeeded = true;
118 mutable MEDCoupling::MCAuto<MEDCoupling::MEDCouplingUMesh> _meshorig;
119 mutable MEDCoupling::MCAuto<MEDCoupling::DataArrayIdType> _untouched_2d_cells;
120 mutable MEDCoupling::MCAuto<MEDCoupling::DataArrayIdType> _cells_at_boundary_origin;
121 mutable MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> _centers;
122 mutable MEDCoupling::MCAuto<MEDCoupling::DataArrayDouble> _vol;
124 std::vector< std::unique_ptr< vtkInternal > > Internal2;
127 vtkSedimentDeposit(const vtkSedimentDeposit &) = delete;
128 void operator=(const vtkSedimentDeposit &) = delete;