/[escript]/trunk/finley/src/CPPAdapter/MeshAdapter.h
ViewVC logotype

Contents of /trunk/finley/src/CPPAdapter/MeshAdapter.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 626 - (show annotations)
Thu Mar 23 02:16:36 2006 UTC (12 years, 9 months ago) by elspeth
File MIME type: text/plain
File size: 11675 byte(s)
Copyright information inserted

1 // $Id$
2 /*
3 ************************************************************
4 * Copyright 2006 by ACcESS MNRF *
5 * *
6 * http://www.access.edu.au *
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 ************************************************************
12 */
13
14 #if !defined finley_MeshAdapter_20040526_H
15 #define finley_MeshAdapter_20040526_H
16
17 extern "C" {
18 #include "Mesh.h"
19 #include "Finley.h"
20 #include "Assemble.h"
21 #include "Finley.h"
22 #include "SystemMatrix.h"
23 }
24
25 #include "FinleyError.h"
26 #include "FinleyAdapterException.h"
27
28 #include "SystemMatrixAdapter.h"
29 #include "AbstractContinuousDomain.h"
30 #include "FunctionSpace.h"
31 #include "FunctionSpaceFactory.h"
32
33 #include <boost/shared_ptr.hpp>
34 #include <boost/python/dict.hpp>
35 #include <boost/python/extract.hpp>
36
37 #include <map>
38 #include <vector>
39 #include <string>
40 #include <sstream>
41
42 //
43 // forward declarations
44 class Data;
45
46 //using namespace escript;
47
48 namespace finley {
49
50 struct null_deleter
51 {
52 void operator()(void const *ptr) const
53 {
54 }
55 };
56
57
58 /**
59 \brief
60 MeshAdapter implements the AbstractContinuousDomain
61 interface for the Finley library.
62
63 Description:
64 MeshAdapter implements the AbstractContinuousDomain
65 interface for the Finley library.
66 */
67
68 class MeshAdapter : public escript::AbstractContinuousDomain {
69
70 public:
71
72 //
73 // Codes for function space types supported
74 static const int DegreesOfFreedom;
75 static const int ReducedDegreesOfFreedom;
76 static const int Nodes;
77 static const int Elements;
78 static const int FaceElements;
79 static const int Points;
80 static const int ContactElementsZero;
81 static const int ContactElementsOne;
82
83 /**
84 \brief
85 Constructor for MeshAdapter
86
87 Description:
88 Constructor for MeshAdapter. The pointer passed to MeshAdapter
89 is deleted using a call to Finley_Mesh_deallocate in the
90 MeshAdapter destructor.
91
92 Throws:
93 May throw an exception derived from EsysException
94
95 \param finleyMesh Input - A pointer to the externally constructed
96 finley mesh.The pointer passed to MeshAdapter
97 is deleted using a call to
98 Finley_Mesh_deallocate in the MeshAdapter
99 destructor.
100 */
101 MeshAdapter(Finley_Mesh* finleyMesh=0);
102
103 /**
104 \brief
105 Copy constructor.
106 */
107 MeshAdapter(const MeshAdapter& in);
108
109 /**
110 \brief
111 Destructor for MeshAdapter. As specified in the constructor
112 this calls Finley_Mesh_deallocate for the pointer given to the
113 constructor.
114 */
115 ~MeshAdapter();
116
117 /**
118 \brief
119 return this as an AbstractContinuousDomain.
120 */
121 inline const AbstractContinuousDomain& asAbstractContinuousDomain() const
122 {
123 return *(static_cast<const AbstractContinuousDomain*>(this));
124 }
125
126 /**
127 \brief
128 return this as an AbstractDomain.
129 */
130 inline const AbstractDomain& asAbstractDomain() const
131 {
132 return *(static_cast<const AbstractDomain*>(this));
133 }
134
135 /**
136 \brief
137 Write the current mesh to a file with the given name.
138 \param fileName Input - The name of the file to write to.
139 */
140 void write(const std::string& fileName) const;
141
142 /**
143 \brief
144 return the pointer to the underlying finley mesh structure
145 */
146 Finley_Mesh* getFinley_Mesh() const;
147
148 /**
149 \brief
150 Return the tag key for the given sample number.
151 \param functionSpaceType Input - The function space type.
152 \param sampleNo Input - The sample number.
153 */
154 int getTagFromSampleNo(int functionSpaceType, int sampleNo) const;
155
156 /**
157 \brief
158 Return the reference number of the given sample number.
159 \param functionSpaceType Input - The function space type.
160 \param sampleNo Input - The sample number.
161 */
162 int getReferenceNoFromSampleNo(int functionSpaceType, int sampleNo) const;
163
164 /**
165 \brief
166 Returns true if the given integer is a valid function space type
167 for this domain.
168 */
169 virtual bool isValidFunctionSpaceType(int functionSpaceType) const;
170
171 /**
172 \brief
173 Return a description for this domain
174 */
175 virtual std::string getDescription() const;
176
177 /**
178 \brief
179 Return a description for the given function space type code
180 */
181 virtual std::string functionSpaceTypeAsString(int functionSpaceType) const;
182
183 /**
184 \brief
185 Build the table of function space type names
186 */
187 void setFunctionSpaceTypeNames();
188
189 /**
190 \brief
191 Return a continuous FunctionSpace code
192 */
193 virtual int getContinuousFunctionCode() const;
194
195 /**
196 \brief
197 Return a functon FunctionSpace code
198 */
199 virtual int getFunctionCode() const;
200
201 /**
202 \brief
203 Return a function on boundary FunctionSpace code
204 */
205 virtual int getFunctionOnBoundaryCode() const;
206
207 /**
208 \brief
209 Return a FunctionOnContactZero code
210 */
211 virtual int getFunctionOnContactZeroCode() const;
212
213 /**
214 \brief
215 Return a FunctionOnContactOne code
216 */
217 virtual int getFunctionOnContactOneCode() const;
218
219 /**
220 \brief
221 Return a Solution code
222 */
223 virtual int getSolutionCode() const;
224
225 /**
226 \brief
227 Return a ReducedSolution code
228 */
229 virtual int getReducedSolutionCode() const;
230
231 /**
232 \brief
233 Return a DiracDeltaFunction code
234 */
235 virtual int getDiracDeltaFunctionCode() const;
236
237 /**
238 \brief
239 */
240 typedef std::map<int, std::string> FunctionSpaceNamesMapType;
241
242 /**
243 \brief
244 */
245 virtual int getDim() const;
246
247 /**
248 \brief
249 Return the number of data points per sample, and the number of samples as a pair.
250 \param functionSpaceCode Input -
251 */
252 virtual std::pair<int,int> getDataShape(int functionSpaceCode) const;
253
254 /**
255 \brief
256 copies the location of data points into arg. The domain of arg has to match this.
257 has to be implemented by the actual Domain adapter.
258 */
259 virtual void setToX(escript::Data& arg) const;
260
261 /**
262 \brief
263 assigns new location to the domain
264 */
265 virtual void setNewX(const escript::Data& arg);
266
267 /**
268 \brief
269 interpolates data given on source onto target where source and target have to be given on the same domain.
270 */
271 virtual void interpolateOnDomain(escript::Data& target,const escript::Data& source) const;
272 virtual bool probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const;
273
274 /**
275 \brief
276 interpolates data given on source onto target where source and target are given on different domains.
277 has to be implemented by the actual Domain adapter.
278 */
279 virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
280 virtual bool probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const;
281
282 /**
283 \brief
284 copies the surface normals at data points into out. The actual function space to be considered
285 is defined by out. out has to be defined on this.
286 */
287 virtual void setToNormal(escript::Data& out) const;
288
289 /**
290 \brief
291 copies the size of samples into out. The actual function space to be considered
292 is defined by out. out has to be defined on this.
293 */
294 virtual void setToSize(escript::Data& out) const;
295
296 /**
297 \brief
298 copies the gradient of arg into grad. The actual function space to be considered
299 for the gradient is defined by grad. arg and grad have to be defined on this.
300 */
301 virtual void setToGradient(escript::Data& grad,const escript::Data& arg) const;
302
303 /**
304 \brief
305 copies the integrals of the function defined by arg into integrals.
306 arg has to be defined on this.
307 */
308 virtual void setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const;
309
310 /**
311 \brief
312 return the identifier of the matrix type to be used for the global stiffness matrix when a particular solver, package
313 and symmetric matrix is used.
314 \param solver
315 \param symmetry
316 */
317 virtual int getSystemMatrixTypeId(const int solver, const int package, const bool symmetry) const;
318
319 /**
320 \brief
321 returns true if data on this domain and a function space of type functionSpaceCode has to
322 considered as cell centered data.
323 */
324 virtual bool isCellOriented(int functionSpaceCode) const;
325
326 /**
327 \brief
328 Saves a dictonary of Data objects to an OpenDX input file. The keywords are used as identifier
329
330 This has to be implemented by the actual Domain adapter.
331 */
332 virtual void saveDX(const std::string& filename,const boost::python::dict& arg) const;
333
334
335 /**
336 \brief
337 Saves a dictonary of Data objects to an VTK XML input file. The keywords are used as identifier
338
339 This has to be implemented by the actual Domain adapter.
340 */
341 virtual void saveVTK(const std::string& filename,const boost::python::dict& arg) const;
342
343 /**
344 \brief
345 returns the function space representation of the type functionSpaceCode on this domain
346 as a vtkObject.
347 */
348 // vtkObject createVtkObject(int functionSpaceCode) const;
349
350 /**
351 \brief
352 adds a PDE onto the stiffness matrix mat and a rhs
353 */
354 virtual void addPDEToSystem(
355 SystemMatrixAdapter& mat, escript::Data& rhs,
356 const escript::Data& A, const escript::Data& B, const escript::Data& C,
357 const escript::Data& D, const escript::Data& X, const escript::Data& Y,
358 const escript::Data& d, const escript::Data& y,
359 const escript::Data& d_contact, const escript::Data& y_contact) const;
360
361 /**
362 \brief
363 adds a PDE onto the stiffness matrix mat and a rhs
364 */
365 virtual void addPDEToRHS(escript::Data& rhs,
366 const escript::Data& X, const escript::Data& Y,
367 const escript::Data& y, const escript::Data& y_contact) const;
368
369 /**
370 \brief
371 creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros:
372 */
373 SystemMatrixAdapter newSystemMatrix(
374 const int row_blocksize,
375 const escript::FunctionSpace& row_functionspace,
376 const int column_blocksize,
377 const escript::FunctionSpace& column_functionspace,
378 const int type) const;
379
380 /**
381 \brief returns locations in the FEM nodes
382 */
383 virtual escript::Data getX() const;
384
385 /**
386 \brief return boundary normals at the quadrature point on the face elements
387 */
388 virtual escript::Data getNormal() const;
389
390 /**
391 \brief returns the element size
392 */
393 virtual escript::Data getSize() const;
394
395 /**
396 \brief comparison operators
397 */
398 virtual bool operator==(const AbstractDomain& other) const;
399 virtual bool operator!=(const AbstractDomain& other) const;
400
401 protected:
402
403 private:
404
405 //
406 // pointer to the externally created finley mesh
407 boost::shared_ptr<Finley_Mesh> m_finleyMesh;
408
409 static FunctionSpaceNamesMapType m_functionSpaceTypeNames;
410
411 };
412
413 } // end of namespace
414
415 #endif

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26