Escript  Revision_4320
Rectangle.h
Go to the documentation of this file.
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_RECTANGLE_H__
17 #define __RIPLEY_RECTANGLE_H__
18 
19 #include <ripley/RipleyDomain.h>
20 
21 struct Paso_Connector;
22 
23 namespace ripley {
24 
30 {
31 public:
32 
40  Rectangle(int n0, int n1, double x0, double y0, double x1, double y1,
41  int d0=-1, int d1=-1);
42 
47  ~Rectangle();
48 
53  virtual std::string getDescription() const;
54 
58  virtual bool operator==(const escript::AbstractDomain& other) const;
59 
65  void dump(const std::string& filename) const;
66 
69  virtual void readBinaryGrid(escript::Data& out, std::string filename,
70  const std::vector<int>& first,
71  const std::vector<int>& numValues,
72  const std::vector<int>& multiplier) const;
73 
76  virtual void readNcGrid(escript::Data& out, std::string filename,
77  std::string varname, const std::vector<int>& first,
78  const std::vector<int>& numValues,
79  const std::vector<int>& multiplier) const;
80 
86  const int* borrowSampleReferenceIDs(int fsType) const;
87 
92  virtual bool ownSample(int fs_code, index_t id) const;
93 
100  virtual void setToNormal(escript::Data& out) const;
101 
107  virtual void setToSize(escript::Data& out) const;
108 
113  virtual int getNumDataPointsGlobal() const { return (m_gNE0+1)*(m_gNE1+1); }
114 
120  virtual void Print_Mesh_Info(const bool full=false) const;
121 
126  virtual IndexVector getNumNodesPerDim() const;
127 
132  virtual IndexVector getNumElementsPerDim() const;
133 
139  virtual IndexVector getNumFacesPerBoundary() const;
140 
145  virtual IndexVector getNodeDistribution() const { return m_nodeDistribution; }
146 
151  virtual IndexVector getNumSubdivisionsPerDim() const;
152 
158  virtual std::pair<double,double> getFirstCoordAndSpacing(dim_t dim) const;
159 
160 protected:
161  virtual dim_t getNumNodes() const { return m_N0*m_N1; }
162  virtual dim_t getNumElements() const { return m_NE0*m_NE1; }
163  virtual dim_t getNumFaceElements() const;
164  virtual dim_t getNumDOF() const;
165  virtual dim_t insertNeighbourNodes(IndexVector& index, index_t node) const;
166  virtual void assembleCoordinates(escript::Data& arg) const;
167  virtual void assembleGradient(escript::Data& out, escript::Data& in) const;
168  virtual void assembleIntegrate(std::vector<double>& integrals, escript::Data& arg) const;
169  virtual void assemblePDESingle(Paso_SystemMatrix* mat, escript::Data& rhs,
170  const escript::Data& A, const escript::Data& B,
171  const escript::Data& C, const escript::Data& D,
172  const escript::Data& X, const escript::Data& Y) const;
173  virtual void assemblePDEBoundarySingle(Paso_SystemMatrix* mat,
174  escript::Data& rhs, const escript::Data& d,
175  const escript::Data& y) const;
176  virtual void assemblePDESingleReduced(Paso_SystemMatrix* mat,
177  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
178  const escript::Data& C, const escript::Data& D,
179  const escript::Data& X, const escript::Data& Y) const;
180  virtual void assemblePDEBoundarySingleReduced(Paso_SystemMatrix* mat,
181  escript::Data& rhs, const escript::Data& d,
182  const escript::Data& y) const;
183  virtual void assemblePDESystem(Paso_SystemMatrix* mat, escript::Data& rhs,
184  const escript::Data& A, const escript::Data& B,
185  const escript::Data& C, const escript::Data& D,
186  const escript::Data& X, const escript::Data& Y) const;
187  virtual void assemblePDEBoundarySystem(Paso_SystemMatrix* mat,
188  escript::Data& rhs, const escript::Data& d,
189  const escript::Data& y) const;
190  virtual void assemblePDESystemReduced(Paso_SystemMatrix* mat,
191  escript::Data& rhs, const escript::Data& A, const escript::Data& B,
192  const escript::Data& C, const escript::Data& D,
193  const escript::Data& X, const escript::Data& Y) const;
194  virtual void assemblePDEBoundarySystemReduced(Paso_SystemMatrix* mat,
195  escript::Data& rhs, const escript::Data& d,
196  const escript::Data& y) const;
197  virtual Paso_SystemMatrixPattern* getPattern(bool reducedRowOrder, bool reducedColOrder) const;
198  virtual void interpolateNodesOnElements(escript::Data& out,
199  escript::Data& in, bool reduced) const;
200  virtual void interpolateNodesOnFaces(escript::Data& out, escript::Data& in,
201  bool reduced) const;
202  virtual void nodesToDOF(escript::Data& out, escript::Data& in) const;
203  virtual void dofToNodes(escript::Data& out, escript::Data& in) const;
204 
205 private:
206  void populateSampleIds();
207  void createPattern();
208  void addToMatrixAndRHS(Paso_SystemMatrix* S, escript::Data& F,
209  const std::vector<double>& EM_S, const std::vector<double>& EM_F,
210  bool addS, bool addF, index_t firstNode, dim_t nEq=1, dim_t nComp=1) const;
211 
213  dim_t m_gNE0, m_gNE1;
214 
216  double m_x0, m_y0;
217 
219  double m_l0, m_l1;
220 
222  int m_NX, m_NY;
223 
225  dim_t m_NE0, m_NE1;
226 
228  dim_t m_ownNE0, m_ownNE1;
229 
231  dim_t m_N0, m_N1;
232 
234  dim_t m_offset0, m_offset1;
235 
239 
245 
246  // vector with first node id on each rank
248 
249  // vector that maps each node to a DOF index (used for the coupler)
251 
252  // Paso connector used by the system matrix and to interpolate DOF to
253  // nodes
255 
256  // the Paso System Matrix pattern
258 };
259 
260 } // end of namespace ripley
261 
262 #endif // __RIPLEY_RECTANGLE_H__
263