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

Contents of /branches/ROBW_XPLATFORM/finley/src/CPPAdapter/MeshAdapter.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 624 - (show annotations)
Wed Mar 22 23:00:01 2006 UTC (15 years, 1 month ago) by robwdcock
File MIME type: text/plain
File size: 12024 byte(s)
+ Finish off the header file path changes
+ The C/C++ libraries now all compile on the ivec altix (cognac) using the xplatform build system

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26