core_example_00001.cpp
1 /*---------------------------------------------------------------------------*\
2  *
3  * mimmo
4  *
5  * Copyright (C) 2015-2021 OPTIMAD engineering Srl
6  *
7  * -------------------------------------------------------------------------
8  * License
9  * This file is part of mimmo.
10  *
11  * mimmo is free software: you can redistribute it and/or modify it
12  * under the terms of the GNU Lesser General Public License v3 (LGPL)
13  * as published by the Free Software Foundation.
14  *
15  * mimmo is distributed in the hope that it will be useful, but WITHOUT
16  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
17  * FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
18  * License for more details.
19  *
20  * You should have received a copy of the GNU Lesser General Public License
21  * along with mimmo. If not, see <http://www.gnu.org/licenses/>.
22  *
23  \ *---------------------------------------------------------------------------*/
24 
25 
26 #include "mimmo_core.hpp"
27 #include "mimmo_iogeneric.hpp"
28 #if MIMMO_ENABLE_MPI
29 #include "mimmo_parallel.hpp"
30 #endif
31 
32 // =================================================================================== //
47 void test00001() {
48 
49  /*
50  Reading a STL geometry from file.
51  */
53  mimmo0->setReadDir("geodata");
54  mimmo0->setReadFileType(FileType::STL);
55  mimmo0->setReadFilename("plane3");
56  mimmo0->execute();
57 
58 #if MIMMO_ENABLE_MPI
59  /*
60  MPI version : distribute the read geometry over X ranks.
61  */
62  mimmo::Partition* partition= new mimmo::Partition();
63  partition->setPartitionMethod(mimmo::PartitionMethod::PARTGEOM);
64  partition->setGeometry(mimmo0->getGeometry());
65  partition->execute();
66 #endif
67 
68  /*
69  * Creation of a synthetic float field on geometry mesh POINT location.
70  */
72  pointField.initialize(mimmo0->getGeometry(), mimmo::MPVLocation::POINT, 1.);
73  std::array<double,3> center({{0.5,0.,0.}});
74  for (bitpit::Vertex & vertex : mimmo0->getGeometry()->getPatch()->getVertices()){
75  std::array<double,3> coords = vertex.getCoords();
76  double value = norm2(coords-center);
77  pointField[vertex.getId()] = value;
78  }
79 
80  /*
81  * Write POINT field and geometry on file.
82  */
83  {
84  std::vector<double> field(mimmo0->getGeometry()->getNVertices());
85  int count = 0;
86  for (bitpit::Vertex & vertex : mimmo0->getGeometry()->getPatch()->getVertices()){
87  field[count] = pointField[vertex.getId()];
88  count++;
89  }
90  mimmo0->getGeometry()->getPatch()->getVTK().addData("pointField", bitpit::VTKFieldType::SCALAR, bitpit::VTKLocation::POINT, field);
91  mimmo0->getGeometry()->getPatch()->write("core_example_00001.0000");
92  mimmo0->getGeometry()->getPatch()->getVTK().removeData("pointField");
93  }
94 
95  /*
96  * Interpolation of synthetic field from POINT to CELL location.
97  */
98  double p = 5.;
100 
101  /*
102  * Write CELL field and geometry on file.
103  */
104  {
105  std::vector<double> field(mimmo0->getGeometry()->getNCells());
106  int count = 0;
107  for (bitpit::Cell & cell : mimmo0->getGeometry()->getPatch()->getCells()){
108  field[count] = cellField[cell.getId()];
109  count++;
110  }
111  mimmo0->getGeometry()->getPatch()->getVTK().addData("cellField", bitpit::VTKFieldType::SCALAR, bitpit::VTKLocation::CELL, field);
112  mimmo0->getGeometry()->getPatch()->write("core_example_00001.0001");
113  mimmo0->getGeometry()->getPatch()->getVTK().removeData("cellField");
114  }
115 
116  /*
117  * Interpolate back CELL synthetic field again on points.
118  */
119  p = 5.;
120  mimmo::MimmoPiercedVector<double> pointField2 = cellField.cellDataToPointData(p);
121 
122  /*
123  * Write new POINT field and geometry on file.
124  */
125  {
126  std::vector<double> field(mimmo0->getGeometry()->getNVertices());
127  int count = 0;
128  for (bitpit::Vertex & vertex : mimmo0->getGeometry()->getPatch()->getVertices()){
129  field[count] = pointField2[vertex.getId()];
130  count++;
131  }
132  mimmo0->getGeometry()->getPatch()->getVTK().addData("pointField", bitpit::VTKFieldType::SCALAR, bitpit::VTKLocation::POINT, field);
133  mimmo0->getGeometry()->getPatch()->write("core_example_00001.0002");
134  mimmo0->getGeometry()->getPatch()->getVTK().removeData("pointField");
135  }
136 
137  /*
138  * Interpolate synthetic POINT field on boundary edges of the geometry mesh.
139  */
140  mimmo0->getGeometry()->updateInterfaces();
141  p = 5.;
143 
144  /*
145  * Write the boundary edge field on screen.
146  */
147  {
148  std::vector<double> field(mimmo0->getGeometry()->getPatch()->getInterfaceCount());
149  int count = 0;
150  for (bitpit::Interface & interface : mimmo0->getGeometry()->getPatch()->getInterfaces()){
151  if (interface.isBorder()){
152  field[count] = interfaceField[interface.getId()];
153  bitpit::log::cout() << bitpit::log::priority(bitpit::log::NORMAL);
154  bitpit::log::cout() << bitpit::log::visibility(bitpit::log::GLOBAL);
155  bitpit::log::cout() << "field on boundary interface " << count << " : " << field[count] << std::endl;
156  count++;
157  }
158  }
159  }
160 
161  /*
162  Clean up & exit;
163  */
164  delete mimmo0;
165 #if MIMMO_ENABLE_MPI
166  delete partition;
167 #endif
168 
169  return;
170 }
171 
172 
173 int main(int argc, char *argv[]) {
174 
175  BITPIT_UNUSED(argc);
176  BITPIT_UNUSED(argv);
177 
178 #if MIMMO_ENABLE_MPI==1
179  MPI_Init(&argc, &argv);
180 
181  {
182 #endif
183 
185  try{
186  test00001() ;
187  }
188  catch(std::exception & e){
189  std::cout<<"core_example_00001 exited with an error of type : "<<e.what()<<std::endl;
190  return 1;
191  }
192 
193 #if MIMMO_ENABLE_MPI==1
194  }
195 
196  MPI_Finalize();
197 #endif
198 
199  return 0;
200 }
MimmoPiercedVector cellDataToPointData(double p=0.)
MimmoPiercedVector pointDataToCellData(double p=0.)
void setReadFilename(std::string filename)
MimmoGeometry is an executable block class wrapping(linking or internally instantiating) a Mimmo Obje...
void initialize(MimmoSharedPointer< MimmoObject >, MPVLocation, const mpv_t &)
void setReadFileType(FileType type)
MimmoSharedPointer< MimmoObject > getGeometry()
void setReadDir(std::string dir)
MimmoPiercedVector pointDataToBoundaryInterfaceData(double p=0.)