/[escript]/trunk/ripley/src/Brick.h
ViewVC logotype

Contents of /trunk/ripley/src/Brick.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6651 - (show annotations)
Wed Feb 7 02:12:08 2018 UTC (19 months, 2 weeks ago) by jfenwick
File MIME type: text/plain
File size: 13187 byte(s)
Make everyone sad by touching all the files

Copyright dates update

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2018 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16
17 #ifndef __RIPLEY_BRICK_H__
18 #define __RIPLEY_BRICK_H__
19
20 #include <ripley/RipleyDomain.h>
21
22 namespace ripley {
23
24 /**
25 \brief
26 Brick is the 3-dimensional implementation of a RipleyDomain.
27 */
28 class RIPLEY_DLL_API Brick: public RipleyDomain
29 {
30 template<class Scalar> friend class DefaultAssembler3D;
31 friend class WaveAssembler3D;
32 friend class LameAssembler3D;
33 public:
34
35 /**
36 \brief creates a hexagonal mesh with n0 x n1 x n2 elements over the
37 brick [x0,x1] x [y0,y1] x [z0,z1].
38 \param n0,n1,n2 number of elements in each dimension
39 \param x0,y0,z0,x1,y1,z1 coordinates of corner nodes of the brick
40 \param d0,d1,d2 number of subdivisions in each dimension
41 */
42 Brick(dim_t n0, dim_t n1, dim_t n2, double x0, double y0, double z0,
43 double x1, double y1, double z1, int d0=-1, int d1=-1, int d2=-1,
44 const std::vector<double>& points = std::vector<double>(),
45 const std::vector<int>& tags = std::vector<int>(),
46 const TagMap& tagnamestonums = TagMap(),
47 escript::SubWorld_ptr w=escript::SubWorld_ptr());
48
49 /**
50 \brief
51 Destructor.
52 */
53 ~Brick();
54
55 /**
56 \brief
57 returns a description for this domain
58 */
59 virtual std::string getDescription() const;
60
61 /**
62 \brief equality operator
63 */
64 virtual bool operator==(const escript::AbstractDomain& other) const;
65
66 /**
67 \brief
68 writes the current mesh to a file with the given name
69 \param filename The name of the file to write to
70 */
71 virtual void write(const std::string& filename) const;
72
73 /**
74 \brief
75 dumps the mesh to a file with the given name
76 \param filename The name of the output file
77 */
78 void dump(const std::string& filename) const;
79
80 /**
81 */
82 virtual void readNcGrid(escript::Data& out, std::string filename,
83 std::string varname, const ReaderParameters& params) const;
84
85 /**
86 */
87 virtual void readBinaryGrid(escript::Data& out, std::string filename,
88 const ReaderParameters& params) const;
89
90 virtual void readBinaryGridFromZipped(escript::Data& out, std::string filename,
91 const ReaderParameters& params) const;
92
93 /**
94 */
95 virtual void writeBinaryGrid(const escript::Data& in,
96 std::string filename,
97 int byteOrder, int dataType) const;
98
99 /**
100 \brief
101 returns the array of reference numbers for a function space type
102 \param fsType The function space type
103 */
104 const dim_t* borrowSampleReferenceIDs(int fsType) const;
105
106 /**
107 \brief
108 returns true if this rank owns the sample id.
109 */
110 virtual bool ownSample(int fsType, index_t id) const;
111
112 /**
113 \brief
114 copies the surface normals at data points into out. The actual function
115 space to be considered is defined by out. out has to be defined on this
116 domain.
117 */
118 virtual void setToNormal(escript::Data& out) const;
119
120 /**
121 \brief
122 copies the size of samples into out. The actual function space to be
123 considered is defined by out. out has to be defined on this domain.
124 */
125 virtual void setToSize(escript::Data& out) const;
126
127 /**
128 \brief
129 returns the number of data points summed across all MPI processes
130 */
131 virtual dim_t getNumDataPointsGlobal() const;
132
133 /**
134 \brief
135 writes information about the mesh to standard output
136 \param full whether to print additional data
137 */
138 virtual void Print_Mesh_Info(const bool full=false) const;
139
140 /**
141 \brief
142 returns the number of nodes per MPI rank in each dimension
143 */
144 virtual const dim_t* getNumNodesPerDim() const { return m_NN; }
145
146 /**
147 \brief
148 returns the number of elements per MPI rank in each dimension
149 */
150 virtual const dim_t* getNumElementsPerDim() const { return m_NE; }
151
152 /**
153 \brief
154 returns the number of face elements in the order
155 (left,right,bottom,top,front,back) on current MPI rank
156 */
157 virtual const dim_t* getNumFacesPerBoundary() const { return m_faceCount; }
158
159 /**
160 \brief
161 returns the node distribution vector
162 */
163 virtual IndexVector getNodeDistribution() const { return m_nodeDistribution; }
164
165 /**
166 \brief
167 returns the number of spatial subdivisions in each dimension
168 */
169 virtual const int* getNumSubdivisionsPerDim() const { return m_NX; }
170
171 /**
172 \brief
173 returns the index'th coordinate value in given dimension for this rank
174 */
175 virtual double getLocalCoordinate(index_t index, int dim) const;
176
177 /**
178 \brief
179 returns the tuple (origin, spacing, number_of_elements)
180 */
181 virtual boost::python::tuple getGridParameters() const;
182
183 /**
184 * \brief
185 Returns a Data object filled with random data passed through filter.
186 */
187 virtual escript::Data randomFill(const escript::DataTypes::ShapeType& shape,
188 const escript::FunctionSpace& what,
189 long seed,
190 const boost::python::tuple& filter) const;
191
192 virtual Assembler_ptr createAssembler(std::string type,
193 const DataMap& options) const;
194 /**
195 \brief
196 returns the lengths of the domain
197 */
198 const double *getLength() const { return m_length; }
199
200 /**
201 \brief
202 returns the lengths of an element
203 */
204 const double *getElementLength() const { return m_dx; }
205
206 /**
207 \brief
208 returns a vector of rank numbers where vec[i]=n means that rank n
209 'owns' element/face element i.
210 */
211 virtual RankVector getOwnerVector(int fsType) const;
212
213 protected:
214 virtual dim_t getNumNodes() const;
215 virtual dim_t getNumElements() const;
216 virtual dim_t getNumFaceElements() const;
217 virtual dim_t getNumDOF() const;
218 virtual dim_t getNumDOFInAxis(unsigned axis) const;
219 virtual dim_t getFirstInDim(unsigned axis) const;
220
221 virtual IndexVector getDiagonalIndices(bool upperOnly) const;
222 virtual void assembleCoordinates(escript::Data& arg) const;
223 virtual void assembleGradient(escript::Data& out,
224 const escript::Data& in) const;
225 virtual void assembleIntegrate(std::vector<real_t>& integrals,
226 const escript::Data& arg) const;
227 virtual void assembleIntegrate(std::vector<cplx_t>& integrals,
228 const escript::Data& arg) const;
229 virtual std::vector<IndexVector> getConnections(bool includeShared=false) const;
230 #ifdef ESYS_HAVE_TRILINOS
231 virtual esys_trilinos::const_TrilinosGraph_ptr getTrilinosGraph() const;
232 #endif
233
234 #ifdef ESYS_HAVE_PASO
235 virtual paso::SystemMatrixPattern_ptr getPasoMatrixPattern(
236 bool reducedRowOrder, bool reducedColOrder) const;
237 #endif
238 virtual void interpolateNodesOnElements(escript::Data& out,
239 const escript::Data& in, bool reduced) const;
240 virtual void interpolateNodesOnFaces(escript::Data& out,
241 const escript::Data& in,
242 bool reduced) const;
243 virtual void nodesToDOF(escript::Data& out, const escript::Data& in) const;
244 virtual dim_t getDofOfNode(dim_t node) const;
245
246 void populateSampleIds();
247 void populateDofMap();
248
249 template<typename Scalar>
250 void assembleGradientImpl(escript::Data& out,
251 const escript::Data& in) const;
252
253 template<typename Scalar>
254 void addToMatrixAndRHS(escript::AbstractSystemMatrix* S, escript::Data& F,
255 const std::vector<Scalar>& EM_S, const std::vector<Scalar>& EM_F,
256 bool addS, bool addF, index_t firstNode, int nEq=1, int nComp=1) const;
257
258 template<typename ValueType>
259 void readBinaryGridImpl(escript::Data& out, const std::string& filename,
260 const ReaderParameters& params) const;
261 #ifdef ESYS_HAVE_BOOST_IO
262 template<typename ValueType>
263 void readBinaryGridZippedImpl(escript::Data& out,
264 const std::string& filename,
265 const ReaderParameters& params) const;
266 #endif
267 template<typename ValueType>
268 void writeBinaryGridImpl(const escript::Data& in,
269 const std::string& filename, int byteOrder) const;
270
271 dim_t findNode(const double *coords) const;
272
273 virtual escript::Data randomFillWorker(
274 const escript::DataTypes::ShapeType& shape, long seed,
275 const boost::python::tuple& filter) const;
276
277 /// total number of elements in each dimension
278 dim_t m_gNE[3];
279
280 /// origin of domain
281 double m_origin[3];
282
283 /// side lengths of domain
284 double m_length[3];
285
286 /// grid spacings / cell sizes of domain
287 double m_dx[3];
288
289 /// number of spatial subdivisions
290 int m_NX[3];
291
292 /// number of elements for this rank in each dimension including shared
293 dim_t m_NE[3];
294
295 /// number of own elements for this rank in each dimension
296 dim_t m_ownNE[3];
297
298 /// number of nodes for this rank in each dimension
299 dim_t m_NN[3];
300
301 /// first node on this rank is at (offset0,offset1,offset2) in global mesh
302 dim_t m_offset[3];
303
304 /// number of face elements per edge (left, right, bottom, top, front, back)
305 dim_t m_faceCount[6];
306
307 /// faceOffset[i]=-1 if face i is not an external face, otherwise it is
308 /// the index of that face (where i: 0=left, 1=right, 2=bottom, 3=top,
309 /// 4=front, 5=back)
310 IndexVector m_faceOffset;
311
312 /// vector of sample reference identifiers
313 IndexVector m_dofId;
314 IndexVector m_nodeId;
315 IndexVector m_elementId;
316 IndexVector m_faceId;
317
318 // vector with first node id on each rank
319 IndexVector m_nodeDistribution;
320
321 // vector that maps each node to a DOF index
322 IndexVector m_dofMap;
323
324 #ifdef ESYS_HAVE_PASO
325 // the Paso System Matrix pattern
326 mutable paso::SystemMatrixPattern_ptr m_pattern;
327 #endif
328
329 #ifdef ESYS_HAVE_TRILINOS
330 /// Trilinos graph structure, cached for efficiency
331 mutable esys_trilinos::const_TrilinosGraph_ptr m_graph;
332 #endif
333 private:
334 template<typename Scalar>
335 void assembleIntegrateImpl(std::vector<Scalar>& integrals, const escript::Data& arg) const;
336
337 template <typename S>
338 void interpolateNodesOnElementsWorker(escript::Data& out,
339 const escript::Data& in, bool reduced, S sentinel) const;
340 template <typename S>
341 void interpolateNodesOnFacesWorker(escript::Data& out,
342 const escript::Data& in,
343 bool reduced, S sentinel) const;
344 };
345
346 ////////////////////////////// inline methods ////////////////////////////////
347 inline dim_t Brick::getDofOfNode(dim_t node) const
348 {
349 return m_dofMap[node];
350 }
351
352 inline dim_t Brick::getNumDataPointsGlobal() const
353 {
354 return (m_gNE[0]+1)*(m_gNE[1]+1)*(m_gNE[2]+1);
355 }
356
357 inline double Brick::getLocalCoordinate(index_t index, int dim) const
358 {
359 ESYS_ASSERT(dim>=0 && dim<3, "'dim' out of bounds");
360 ESYS_ASSERT(index>=0 && index<m_NN[dim], "'index' out of bounds");
361 return m_origin[dim]+m_dx[dim]*(m_offset[dim]+index);
362 }
363
364 inline boost::python::tuple Brick::getGridParameters() const
365 {
366 return boost::python::make_tuple(
367 boost::python::make_tuple(m_origin[0], m_origin[1], m_origin[2]),
368 boost::python::make_tuple(m_dx[0], m_dx[1], m_dx[2]),
369 boost::python::make_tuple(m_gNE[0], m_gNE[1], m_gNE[2]));
370 }
371
372 //protected
373 inline dim_t Brick::getNumDOF() const
374 {
375 return (m_gNE[0]+1)/m_NX[0]*(m_gNE[1]+1)/m_NX[1]*(m_gNE[2]+1)/m_NX[2];
376 }
377
378 //protected
379 inline dim_t Brick::getNumDOFInAxis(unsigned axis) const
380 {
381 ESYS_ASSERT(axis < m_numDim, "Invalid axis");
382 return (m_gNE[axis]+1)/m_NX[axis];
383 }
384
385 //protected
386 inline dim_t Brick::getNumNodes() const
387 {
388 return m_NN[0]*m_NN[1]*m_NN[2];
389 }
390
391 //protected
392 inline dim_t Brick::getNumElements() const
393 {
394 return m_NE[0]*m_NE[1]*m_NE[2];
395 }
396
397 //protected
398 inline dim_t Brick::getNumFaceElements() const
399 {
400 return m_faceCount[0] + m_faceCount[1] + m_faceCount[2]
401 + m_faceCount[3] + m_faceCount[4] + m_faceCount[5];
402 }
403
404 //protected
405 inline index_t Brick::getFirstInDim(unsigned axis) const
406 {
407 return m_offset[axis] == 0 ? 0 : 1;
408 }
409
410 } // end of namespace ripley
411
412 #endif // __RIPLEY_BRICK_H__
413

  ViewVC Help
Powered by ViewVC 1.1.26