Nektar++
ProcessInterpPointDataToFld.cpp
Go to the documentation of this file.
1 ////////////////////////////////////////////////////////////////////////////////
2 //
3 // File: ProcessInterpPointDataToFld.cpp
4 //
5 // For more information, please see: http://www.nektar.info/
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // Permission is hereby granted, free of charge, to any person obtaining a
14 // copy of this software and associated documentation files (the "Software"),
15 // to deal in the Software without restriction, including without limitation
16 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
17 // and/or sell copies of the Software, and to permit persons to whom the
18 // Software is furnished to do so, subject to the following conditions:
19 //
20 // The above copyright notice and this permission notice shall be included
21 // in all copies or substantial portions of the Software.
22 //
23 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29 // DEALINGS IN THE SOFTWARE.
30 //
31 // Description: Interpolate, using finite different approximation,
32 // given data to a fld file
33 //
34 ////////////////////////////////////////////////////////////////////////////////
35 #include <iostream>
36 #include <string>
37 using namespace std;
38 
39 #include <boost/core/ignore_unused.hpp>
40 #include <boost/geometry.hpp>
41 #include <boost/math/special_functions/fpclassify.hpp>
42 
48 
50 
51 namespace bg = boost::geometry;
52 namespace bgi = boost::geometry::index;
53 
54 namespace Nektar
55 {
56 namespace FieldUtils
57 {
58 
59 ModuleKey ProcessInterpPointDataToFld::className =
61  ModuleKey(eProcessModule, "interppointdatatofld"),
62  ProcessInterpPointDataToFld::create,
63  "Interpolates given discrete data using a finite difference "
64  "approximation to a fld file given a xml file");
65 
66 ProcessInterpPointDataToFld::ProcessInterpPointDataToFld(FieldSharedPtr f)
67  : ProcessModule(f)
68 {
69  m_config["frompts"] = ConfigOption(
70  false, "NotSet", "Pts file from which to interpolate field");
71 
72  m_config["interpcoord"] =
73  ConfigOption(false, "-1", "coordinate id to use for interpolation");
74 }
75 
77 {
78 }
79 
80 void ProcessInterpPointDataToFld::Process(po::variables_map &vm)
81 {
82  m_f->SetUpExp(vm);
83 
84  int i, j;
86  // Load pts file
87  ASSERTL0(m_config["frompts"].as<string>().compare("NotSet") != 0,
88  "ProcessInterpPointDataToFld requires frompts parameter");
89  string inFile = m_config["frompts"].as<string>().c_str();
90 
91  int totpoints = m_f->m_exp[0]->GetTotPoints();
92 
94  for (int i = 0; i < 3; ++i)
95  {
96  intFields[i] = Array<OneD, NekDouble>(totpoints, 0.);
97  }
98  m_f->m_exp[0]->GetCoords(intFields[0], intFields[1], intFields[2]);
99 
100  if (boost::filesystem::path(inFile).extension() == ".pts")
101  {
104 
105  ptsIO->Import(inFile, fieldPts);
106  }
107  else if (boost::filesystem::path(inFile).extension() == ".csv")
108  {
111 
114 
115  NekDouble vmax, vmin, margin = 0.1;
116  vmax = intFields[0][Vmath::Imax(totpoints, intFields[0], 1)];
117  vmin = intFields[0][Vmath::Imin(totpoints, intFields[0], 1)];
118  Range->m_xmax = (vmax - vmin) * margin + vmax;
119  Range->m_xmin = -(vmax - vmin) * margin + vmin;
120 
121  vmax = intFields[1][Vmath::Imax(totpoints, intFields[1], 1)];
122  vmin = intFields[1][Vmath::Imin(totpoints, intFields[1], 1)];
123  Range->m_ymax = (vmax - vmin) * margin + vmax;
124  Range->m_ymin = -(vmax - vmin) * margin + vmin;
125 
126  vmax = intFields[2][Vmath::Imax(totpoints, intFields[2], 1)];
127  vmin = intFields[2][Vmath::Imin(totpoints, intFields[2], 1)];
128  Range->m_zmax = (vmax - vmin) * margin + vmax;
129  Range->m_zmin = -(vmax - vmin) * margin + vmin;
130 
131  csvIO->Import(inFile, fieldPts, LibUtilities::NullFieldMetaDataMap,
132  Range);
133  }
134  else
135  {
136  ASSERTL0(false, "unknown frompts file type");
137  }
138 
139  int nFields = fieldPts->GetNFields();
140  ASSERTL0(nFields > 0, "No field values provided in input");
141 
142  // Define new expansions.
143  ASSERTL0(
144  m_f->m_numHomogeneousDir == 0,
145  "ProcessInterpPointDataToFld does not support homogeneous expansion");
146 
147  m_f->m_exp.resize(nFields);
148  for (i = 1; i < nFields; ++i)
149  {
150  m_f->m_exp[i] = m_f->AppendExpList(m_f->m_numHomogeneousDir);
151  }
152 
153  Array<OneD, Array<OneD, NekDouble>> intFields1(3 + nFields);
154 
155  for (int i = 0; i < 3; ++i)
156  {
157  intFields1[i] = intFields[i];
158  }
159 
160  // Declare space for interpolated fields
161  for (int i = 3; i < 3 + nFields; ++i)
162  {
163  intFields1[i] = Array<OneD, NekDouble>(totpoints);
164  }
165 
168 
169  int coord_id = m_config["interpcoord"].as<int>();
170  ASSERTL0(coord_id <= static_cast<int>(outPts->GetDim()) - 1,
171  "interpcoord is bigger than the Pts files dimension");
172 
173  Interpolator interp(LibUtilities::eNoMethod, coord_id);
174 
175  if (m_f->m_verbose && m_f->m_comm->TreatAsRankZero())
176  {
177  interp.SetProgressCallback(
179  }
180  interp.Interpolate(fieldPts, outPts);
181  if (m_f->m_verbose && m_f->m_comm->TreatAsRankZero())
182  {
183  cout << endl;
184  }
185 
186  for (i = 0; i < totpoints; ++i)
187  {
188  for (j = 0; j < nFields; ++j)
189  {
190  m_f->m_exp[j]->SetPhys(i, outPts->GetPointVal(3 + j, i));
191  }
192  }
193 
194  // forward transform fields
195  for (i = 0; i < nFields; ++i)
196  {
197  m_f->m_exp[i]->FwdTransLocalElmt(m_f->m_exp[i]->GetPhys(),
198  m_f->m_exp[i]->UpdateCoeffs());
199  }
200 
201  // save field names
202  for (int j = 0; j < fieldPts->GetNFields(); ++j)
203  {
204  m_f->m_variables.push_back(fieldPts->GetFieldName(j));
205  }
206 }
207 } // namespace FieldUtils
208 } // namespace Nektar
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215
A class that contains algorithms for interpolation between pts fields, expansions and different meshe...
FIELD_UTILS_EXPORT void Interpolate(const std::vector< MultiRegions::ExpListSharedPtr > expInField, std::vector< MultiRegions::ExpListSharedPtr > &expOutField, NekDouble def_value=0.0)
Interpolate from an expansion to an expansion.
FieldSharedPtr m_f
Field object.
Definition: Module.h:225
std::map< std::string, ConfigOption > m_config
List of configuration values.
Definition: Module.h:228
virtual void Process(po::variables_map &vm)
Write mesh to output file.
void PrintProgressbar(const int position, const int goal) const
Abstract base class for processing modules.
Definition: Module.h:260
void SetProgressCallback(FuncPointerT func, ObjectPointerT obj)
sets a callback funtion which gets called every time the interpolation progresses
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:198
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< Field > FieldSharedPtr
Definition: Field.hpp:989
std::pair< ModuleType, std::string > ModuleKey
Definition: Module.h:285
ModuleFactory & GetModuleFactory()
Definition: Module.cpp:49
static FieldMetaDataMap NullFieldMetaDataMap
Definition: FieldIO.h:53
std::shared_ptr< PtsField > PtsFieldSharedPtr
Definition: PtsField.h:190
std::shared_ptr< DomainRange > DomainRangeShPtr
Definition: DomainRange.h:66
std::shared_ptr< CsvIO > CsvIOSharedPtr
Definition: CsvIO.h:83
std::shared_ptr< PtsIO > PtsIOSharedPtr
Definition: PtsIO.h:97
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:1
double NekDouble
int Imax(int n, const T *x, const int incx)
Return the index of the maximum element in x.
Definition: Vmath.cpp:918
int Imin(int n, const T *x, const int incx)
Return the index of the minimum element in x.
Definition: Vmath.cpp:1023
Represents a command-line configuration option.
Definition: Module.h:131