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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4154 - (show annotations)
Tue Jan 22 09:30:23 2013 UTC (6 years, 7 months ago) by jfenwick
File MIME type: text/plain
File size: 9093 byte(s)
Round 1 of copyright fixes
1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2013 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development since 2012 by School of Earth Sciences
13 *
14 *****************************************************************************/
15
16 #ifndef __RIPLEY_BRICK_H__
17 #define __RIPLEY_BRICK_H__
18
19 #include <ripley/RipleyDomain.h>
20
21 struct Paso_Connector;
22
23 namespace ripley {
24
25 /**
26 \brief
27 Brick is the 3-dimensional implementation of a RipleyDomain.
28 */
29 class RIPLEY_DLL_API Brick: public RipleyDomain
30 {
31 public:
32
33 /**
34 \brief creates a hexagonal mesh with n0 x n1 x n2 elements over the
35 brick [x0,x1] x [y0,y1] x [z0,z1].
36 \param n0,n1,n2 number of elements in each dimension
37 \param x0,y0,z0,x1,y1,z1 coordinates of corner nodes of the brick
38 \param d0,d1,d2 number of subdivisions in each dimension
39 */
40 Brick(int n0, int n1, int n2, double x0, double y0, double z0, double x1,
41 double y1, double z1, int d0=-1, int d1=-1, int d2=-1);
42
43 /**
44 \brief
45 Destructor.
46 */
47 ~Brick();
48
49 /**
50 \brief
51 returns a description for this domain
52 */
53 virtual std::string getDescription() const;
54
55 /**
56 \brief equality operator
57 */
58 virtual bool operator==(const escript::AbstractDomain& other) const;
59
60 /**
61 \brief
62 dumps the mesh to a file with the given name
63 \param filename The name of the output file
64 */
65 void dump(const std::string& filename) const;
66
67 /**
68 */
69 virtual void readBinaryGrid(escript::Data& out, std::string filename,
70 const std::vector<int>& first,
71 const std::vector<int>& numValues) const;
72
73 /**
74 */
75 virtual void readNcGrid(escript::Data& out, std::string filename,
76 std::string varname, const std::vector<int>& first,
77 const std::vector<int>& numValues) const;
78
79 /**
80 \brief
81 returns the reference number of the given sample number
82 \param fsType The function space type
83 */
84 const int* borrowSampleReferenceIDs(int fsType) const;
85
86 /**
87 \brief
88 returns true if this rank owns the sample id.
89 */
90 virtual bool ownSample(int fsType, index_t id) const;
91
92 /**
93 \brief
94 copies the surface normals at data points into out. The actual function
95 space to be considered is defined by out. out has to be defined on this
96 domain.
97 */
98 virtual void setToNormal(escript::Data& out) const;
99
100 /**
101 \brief
102 copies the size of samples into out. The actual function space to be
103 considered is defined by out. out has to be defined on this domain.
104 */
105 virtual void setToSize(escript::Data& out) const;
106
107 /**
108 \brief
109 returns the number of data points summed across all MPI processes
110 */
111 virtual int getNumDataPointsGlobal() const { return (m_gNE0+1)*(m_gNE1+1)*(m_gNE2+1); }
112
113 /**
114 \brief
115 writes information about the mesh to standard output
116 \param full whether to print additional data
117 */
118 virtual void Print_Mesh_Info(const bool full=false) const;
119
120 /**
121 \brief
122 returns the number of nodes per MPI rank in each dimension
123 */
124 virtual IndexVector getNumNodesPerDim() const;
125
126 /**
127 \brief
128 returns the number of elements per MPI rank in each dimension
129 */
130 virtual IndexVector getNumElementsPerDim() const;
131
132 /**
133 \brief
134 returns the number of face elements in the order
135 (left,right,bottom,top,[front,back]) on current MPI rank
136 */
137 virtual IndexVector getNumFacesPerBoundary() const;
138
139 /**
140 \brief
141 returns the node distribution vector
142 */
143 virtual IndexVector getNodeDistribution() const { return m_nodeDistribution; }
144
145 /**
146 \brief
147 returns the number of spatial subdivisions in each dimension
148 */
149 virtual IndexVector getNumSubdivisionsPerDim() const;
150
151 /**
152 \brief
153 returns the first coordinate value and the node spacing along given
154 dimension as a pair
155 */
156 virtual std::pair<double,double> getFirstCoordAndSpacing(dim_t dim) const;
157
158 protected:
159 virtual dim_t getNumNodes() const { return m_N0*m_N1*m_N2; }
160 virtual dim_t getNumElements() const { return m_NE0*m_NE1*m_NE2; }
161 virtual dim_t getNumFaceElements() const;
162 virtual dim_t getNumDOF() const;
163 virtual dim_t insertNeighbourNodes(IndexVector& index, index_t node) const;
164 virtual void assembleCoordinates(escript::Data& arg) const;
165 virtual void assembleGradient(escript::Data& out, escript::Data& in) const;
166 virtual void assembleIntegrate(std::vector<double>& integrals, escript::Data& arg) const;
167 virtual void assemblePDESingle(Paso_SystemMatrix* mat, escript::Data& rhs,
168 const escript::Data& A, const escript::Data& B,
169 const escript::Data& C, const escript::Data& D,
170 const escript::Data& X, const escript::Data& Y) const;
171 virtual void assemblePDEBoundarySingle(Paso_SystemMatrix* mat,
172 escript::Data& rhs, const escript::Data& d,
173 const escript::Data& y) const;
174 virtual void assemblePDESingleReduced(Paso_SystemMatrix* mat,
175 escript::Data& rhs, const escript::Data& A, const escript::Data& B,
176 const escript::Data& C, const escript::Data& D,
177 const escript::Data& X, const escript::Data& Y) const;
178 virtual void assemblePDEBoundarySingleReduced(Paso_SystemMatrix* mat,
179 escript::Data& rhs, const escript::Data& d,
180 const escript::Data& y) const;
181 virtual void assemblePDESystem(Paso_SystemMatrix* mat, escript::Data& rhs,
182 const escript::Data& A, const escript::Data& B,
183 const escript::Data& C, const escript::Data& D,
184 const escript::Data& X, const escript::Data& Y) const;
185 virtual void assemblePDEBoundarySystem(Paso_SystemMatrix* mat,
186 escript::Data& rhs, const escript::Data& d,
187 const escript::Data& y) const;
188 virtual void assemblePDESystemReduced(Paso_SystemMatrix* mat,
189 escript::Data& rhs, const escript::Data& A, const escript::Data& B,
190 const escript::Data& C, const escript::Data& D,
191 const escript::Data& X, const escript::Data& Y) const;
192 virtual void assemblePDEBoundarySystemReduced(Paso_SystemMatrix* mat,
193 escript::Data& rhs, const escript::Data& d,
194 const escript::Data& y) const;
195 virtual Paso_SystemMatrixPattern* getPattern(bool reducedRowOrder, bool reducedColOrder) const;
196 virtual void interpolateNodesOnElements(escript::Data& out,
197 escript::Data& in, bool reduced) const;
198 virtual void interpolateNodesOnFaces(escript::Data& out, escript::Data& in,
199 bool reduced) const;
200 virtual void nodesToDOF(escript::Data& out, escript::Data& in) const;
201 virtual void dofToNodes(escript::Data& out, escript::Data& in) const;
202
203 private:
204 void populateSampleIds();
205 void createPattern();
206 void addToMatrixAndRHS(Paso_SystemMatrix* S, escript::Data& F,
207 const std::vector<double>& EM_S, const std::vector<double>& EM_F,
208 bool addS, bool addF, index_t firstNode, dim_t nEq=1, dim_t nComp=1) const;
209
210
211 /// total number of elements in each dimension
212 dim_t m_gNE0, m_gNE1, m_gNE2;
213
214 /// location of domain
215 double m_x0, m_y0, m_z0;
216
217 /// side lengths of domain
218 double m_l0, m_l1, m_l2;
219
220 /// number of spatial subdivisions
221 int m_NX, m_NY, m_NZ;
222
223 /// number of elements for this rank in each dimension including shared
224 dim_t m_NE0, m_NE1, m_NE2;
225
226 /// number of own elements for this rank in each dimension
227 dim_t m_ownNE0, m_ownNE1, m_ownNE2;
228
229 /// number of nodes for this rank in each dimension
230 dim_t m_N0, m_N1, m_N2;
231
232 /// first node on this rank is at (offset0,offset1,offset2) in global mesh
233 dim_t m_offset0, m_offset1, m_offset2;
234
235 /// faceOffset[i]=-1 if face i is not an external face, otherwise it is
236 /// the index of that face (where i: 0=left, 1=right, 2=bottom, 3=top,
237 /// 4=front, 5=back)
238 IndexVector m_faceOffset;
239
240 /// vector of sample reference identifiers
241 IndexVector m_dofId;
242 IndexVector m_nodeId;
243 IndexVector m_elementId;
244 IndexVector m_faceId;
245
246 // vector with first node id on each rank
247 IndexVector m_nodeDistribution;
248
249 // vector that maps each node to a DOF index (used for the coupler)
250 IndexVector m_dofMap;
251
252 // Paso connector used by the system matrix and to interpolate DOF to
253 // nodes
254 Paso_Connector* m_connector;
255
256 // the Paso System Matrix pattern
257 Paso_SystemMatrixPattern* m_pattern;
258 };
259
260 } // end of namespace ripley
261
262 #endif // __RIPLEY_BRICK_H__
263

  ViewVC Help
Powered by ViewVC 1.1.26