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

Contents of /trunk-mpi-branch/finley/src/CPPAdapter/MeshAdapter.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1223 - (show annotations)
Fri Aug 3 02:40:39 2007 UTC (11 years, 8 months ago) by gross
File MIME type: text/plain
File size: 14976 byte(s)
first attemt towards an improved MPI version.  

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 #include "system_dep.h"
17
18 extern "C" {
19 #include "../Mesh.h"
20 #include "../Finley.h"
21 #include "../Assemble.h"
22 #include "paso/SystemMatrix.h"
23 }
24
25 #include "FinleyError.h"
26 #include "FinleyAdapterException.h"
27
28 #include "SystemMatrixAdapter.h"
29 #include "escript/AbstractContinuousDomain.h"
30 #include "escript/FunctionSpace.h"
31 #include "escript/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 ReducedNodes;
78 static const int Elements;
79 static const int ReducedElements;
80 static const int FaceElements;
81 static const int ReducedFaceElements;
82 static const int Points;
83 static const int ContactElementsZero;
84 static const int ReducedContactElementsZero;
85 static const int ContactElementsOne;
86 static const int ReducedContactElementsOne;
87
88 /**
89 \brief
90 Constructor for MeshAdapter
91
92 Description:
93 Constructor for MeshAdapter. The pointer passed to MeshAdapter
94 is deleted using a call to Finley_Mesh_deallocate in the
95 MeshAdapter destructor.
96
97 Throws:
98 May throw an exception derived from EsysException
99
100 \param finleyMesh Input - A pointer to the externally constructed
101 finley mesh.The pointer passed to MeshAdapter
102 is deleted using a call to
103 Finley_Mesh_deallocate in the MeshAdapter
104 destructor.
105 */
106 FINLEY_DLL_API
107 MeshAdapter(Finley_Mesh* finleyMesh=0);
108
109 /**
110 \brief
111 Copy constructor.
112 */
113 FINLEY_DLL_API
114 MeshAdapter(const MeshAdapter& in);
115
116 /**
117 \brief
118 Destructor for MeshAdapter. As specified in the constructor
119 this calls Finley_Mesh_deallocate for the pointer given to the
120 constructor.
121 */
122 FINLEY_DLL_API
123 ~MeshAdapter();
124
125 /**
126 \brief
127 return the number of processors used for this domain
128 */
129 FINLEY_DLL_API
130 virtual int getMPISize() const;
131 /**
132 \brief
133 return the number MPI rank of this processor
134 */
135
136 FINLEY_DLL_API
137 virtual int getMPIRank() const;
138
139 /**
140 \brief
141 return this as an AbstractContinuousDomain.
142 */
143 inline const AbstractContinuousDomain& asAbstractContinuousDomain() const
144 {
145 return *(static_cast<const AbstractContinuousDomain*>(this));
146 }
147
148 /**
149 \brief
150 return this as an AbstractDomain.
151 */
152 inline const AbstractDomain& asAbstractDomain() const
153 {
154 return *(static_cast<const AbstractDomain*>(this));
155 }
156
157 /**
158 \brief
159 Write the current mesh to a file with the given name.
160 \param fileName Input - The name of the file to write to.
161 */
162 FINLEY_DLL_API
163 void write(const std::string& fileName) const;
164
165 /**
166 \brief
167 return the pointer to the underlying finley mesh structure
168 */
169 FINLEY_DLL_API
170 Finley_Mesh* getFinley_Mesh() const;
171
172 /**
173 \brief
174 Return the tag key for the given sample number.
175 \param functionSpaceType Input - The function space type.
176 \param sampleNo Input - The sample number.
177 */
178 FINLEY_DLL_API
179 int getTagFromSampleNo(int functionSpaceType, int sampleNo) const;
180
181 /**
182 \brief
183 Return the reference number of the given sample number.
184 \param functionSpaceType Input - The function space type.
185 */
186 FINLEY_DLL_API
187 int* borrowSampleReferenceIDs(int functionSpaceType) const;
188
189 /**
190 \brief
191 Returns true if the given integer is a valid function space type
192 for this domain.
193 */
194 FINLEY_DLL_API
195 virtual bool isValidFunctionSpaceType(int functionSpaceType) const;
196
197 /**
198 \brief
199 Return a description for this domain
200 */
201 FINLEY_DLL_API
202 virtual std::string getDescription() const;
203
204 /**
205 \brief
206 Return a description for the given function space type code
207 */
208 FINLEY_DLL_API
209 virtual std::string functionSpaceTypeAsString(int functionSpaceType) const;
210
211 /**
212 \brief
213 Build the table of function space type names
214 */
215 FINLEY_DLL_API
216 void setFunctionSpaceTypeNames();
217
218 /**
219 \brief
220 Return a continuous FunctionSpace code
221 */
222 FINLEY_DLL_API
223 virtual int getContinuousFunctionCode() const;
224
225 /**
226 \brief
227 Return a continuous on reduced order nodes FunctionSpace code
228 */
229 FINLEY_DLL_API
230 virtual int getReducedContinuousFunctionCode() const;
231
232 /**
233 \brief
234 Return a function FunctionSpace code
235 */
236 FINLEY_DLL_API
237 virtual int getFunctionCode() const;
238
239 /**
240 \brief
241 Return a function with reduced integration order FunctionSpace code
242 */
243 FINLEY_DLL_API
244 virtual int getReducedFunctionCode() const;
245
246 /**
247 \brief
248 Return a function on boundary FunctionSpace code
249 */
250 FINLEY_DLL_API
251 virtual int getFunctionOnBoundaryCode() const;
252
253 /**
254 \brief
255 Return a function on boundary with reduced integration order FunctionSpace code
256 */
257 FINLEY_DLL_API
258 virtual int getReducedFunctionOnBoundaryCode() const;
259
260 /**
261 \brief
262 Return a FunctionOnContactZero code
263 */
264 FINLEY_DLL_API
265 virtual int getFunctionOnContactZeroCode() const;
266
267 /**
268 \brief
269 Return a FunctionOnContactZero code with reduced integration order
270 */
271 FINLEY_DLL_API
272 virtual int getReducedFunctionOnContactZeroCode() const;
273
274 /**
275 \brief
276 Return a FunctionOnContactOne code
277 */
278 FINLEY_DLL_API
279 virtual int getFunctionOnContactOneCode() const;
280
281 /**
282 \brief
283 Return a FunctionOnContactOne code with reduced integration order
284 */
285 FINLEY_DLL_API
286 virtual int getReducedFunctionOnContactOneCode() const;
287
288 /**
289 \brief
290 Return a Solution code
291 */
292 FINLEY_DLL_API
293 virtual int getSolutionCode() const;
294
295 /**
296 \brief
297 Return a ReducedSolution code
298 */
299 FINLEY_DLL_API
300 virtual int getReducedSolutionCode() const;
301
302 /**
303 \brief
304 Return a DiracDeltaFunction code
305 */
306 FINLEY_DLL_API
307 virtual int getDiracDeltaFunctionCode() const;
308
309 /**
310 5B
311 \brief
312 */
313 typedef std::map<int, std::string> FunctionSpaceNamesMapType;
314
315 /**
316 \brief
317 */
318 FINLEY_DLL_API
319 virtual int getDim() const;
320
321 /**
322 \brief
323 Return the number of data points per sample, and the number of samples as a pair.
324 \param functionSpaceCode Input -
325 */
326 FINLEY_DLL_API
327 virtual std::pair<int,int> getDataShape(int functionSpaceCode) const;
328
329 /**
330 \brief
331 copies the location of data points into arg. The domain of arg has to match this.
332 has to be implemented by the actual Domain adapter.
333 */
334 FINLEY_DLL_API
335 virtual void setToX(escript::Data& arg) const;
336
337 /**
338 \brief
339 sets a map from a clear tag name to a tag key
340 \param name Input - tag name.
341 \param tag Input - tag key.
342 */
343 FINLEY_DLL_API
344 virtual void setTagMap(const std::string& name, int tag);
345
346 /**
347 \brief
348 Return the tag key for tag name.
349 \param name Input - tag name
350 */
351 FINLEY_DLL_API
352 virtual int getTag(const std::string& name) const;
353
354 /**
355 \brief
356 Returns true if name is a defined tage name.
357 \param name Input - tag name to be checked.
358 */
359 FINLEY_DLL_API
360 virtual bool isValidTagName(const std::string& name) const;
361
362 /**
363 \brief
364 Returns all tag names in a single string sperated by commas
365 */
366 FINLEY_DLL_API
367 virtual std::string showTagNames() const;
368
369 /**
370 \brief
371 assigns new location to the domain
372 */
373 FINLEY_DLL_API
374 virtual void setNewX(const escript::Data& arg);
375
376 /**
377 \brief
378 interpolates data given on source onto target where source and target have to be given on the same domain.
379 */
380 FINLEY_DLL_API
381 virtual void interpolateOnDomain(escript::Data& target,const escript::Data& source) const;
382 FINLEY_DLL_API
383 virtual bool probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const;
384
385 /**
386 \brief
387 interpolates data given on source onto target where source and target are given on different domains.
388 has to be implemented by the actual Domain adapter.
389 */
390 FINLEY_DLL_API
391 virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
392 FINLEY_DLL_API
393 virtual bool probeInterpolationACross(int functionSpaceType_source,const AbstractDomain& targetDomain, int functionSpaceType_target) const;
394
395 /**
396 \brief
397 copies the surface normals at data points into out. The actual function space to be considered
398 is defined by out. out has to be defined on this.
399 */
400 FINLEY_DLL_API
401 virtual void setToNormal(escript::Data& out) const;
402
403 /**
404 \brief
405 copies the size of samples into out. The actual function space to be considered
406 is defined by out. out has to be defined on this.
407 */
408 FINLEY_DLL_API
409 virtual void setToSize(escript::Data& out) const;
410
411 /**
412 \brief
413 copies the gradient of arg into grad. The actual function space to be considered
414 for the gradient is defined by grad. arg and grad have to be defined on this.
415 */
416 FINLEY_DLL_API
417 virtual void setToGradient(escript::Data& grad,const escript::Data& arg) const;
418
419 /**
420 \brief
421 copies the integrals of the function defined by arg into integrals.
422 arg has to be defined on this.
423 */
424 FINLEY_DLL_API
425 virtual void setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const;
426
427 /**
428 \brief
429 return the identifier of the matrix type to be used for the global stiffness matrix when a particular solver, package
430 and symmetric matrix is used.
431 \param solver
432 \param symmetry
433 */
434 FINLEY_DLL_API
435 virtual int getSystemMatrixTypeId(const int solver, const int package, const bool symmetry) const;
436
437 /**
438 \brief
439 returns true if data on this domain and a function space of type functionSpaceCode has to
440 considered as cell centered data.
441 */
442 FINLEY_DLL_API
443 virtual bool isCellOriented(int functionSpaceCode) const;
444
445 /**
446 \brief
447 Saves a dictonary of Data objects to an OpenDX input file. The keywords are used as identifier
448
449 This has to be implemented by the actual Domain adapter.
450 */
451 FINLEY_DLL_API
452 virtual void saveDX(const std::string& filename,const boost::python::dict& arg) const;
453
454
455 /**
456 \brief
457 Saves a dictonary of Data objects to an VTK XML input file. The keywords are used as identifier
458
459 This has to be implemented by the actual Domain adapter.
460 */
461 FINLEY_DLL_API
462 virtual void saveVTK(const std::string& filename,const boost::python::dict& arg) const;
463
464 /**
465 \brief
466 returns the function space representation of the type functionSpaceCode on this domain
467 as a vtkObject.
468 */
469 // vtkObject createVtkObject(int functionSpaceCode) const;
470
471 /**
472 \brief
473 adds a PDE onto the stiffness matrix mat and a rhs
474 */
475 FINLEY_DLL_API
476 virtual void addPDEToSystem(
477 SystemMatrixAdapter& mat, escript::Data& rhs,
478 const escript::Data& A, const escript::Data& B, const escript::Data& C,
479 const escript::Data& D, const escript::Data& X, const escript::Data& Y,
480 const escript::Data& d, const escript::Data& y,
481 const escript::Data& d_contact, const escript::Data& y_contact) const;
482 /**
483 \brief
484 adds a PDE onto the lumped stiffness matrix matrix
485 */
486 FINLEY_DLL_API
487 virtual void addPDEToLumpedSystem(
488 escript::Data& mat,
489 const escript::Data& D,
490 const escript::Data& d) const;
491
492 /**
493 \brief
494 adds a PDE onto the stiffness matrix mat and a rhs
495 */
496 FINLEY_DLL_API
497 virtual void addPDEToRHS(escript::Data& rhs,
498 const escript::Data& X, const escript::Data& Y,
499 const escript::Data& y, const escript::Data& y_contact) const;
500
501 /**
502 \brief
503 creates a SystemMatrixAdapter stiffness matrix an initializes it with zeros:
504 */
505 FINLEY_DLL_API
506 SystemMatrixAdapter newSystemMatrix(
507 const int row_blocksize,
508 const escript::FunctionSpace& row_functionspace,
509 const int column_blocksize,
510 const escript::FunctionSpace& column_functionspace,
511 const int type) const;
512
513 /**
514 \brief returns locations in the FEM nodes
515 */
516 FINLEY_DLL_API
517 virtual escript::Data getX() const;
518
519 /**
520 \brief return boundary normals at the quadrature point on the face elements
521 */
522 FINLEY_DLL_API
523 virtual escript::Data getNormal() const;
524
525 /**
526 \brief returns the element size
527 */
528 FINLEY_DLL_API
529 virtual escript::Data getSize() const;
530
531 /**
532 \brief comparison operators
533 */
534 FINLEY_DLL_API
535 virtual bool operator==(const AbstractDomain& other) const;
536 FINLEY_DLL_API
537 virtual bool operator!=(const AbstractDomain& other) const;
538
539 /**
540 \brief assigns new tag newTag to all samples of functionspace with a positive
541 value of mask for any its sample point.
542
543 */
544 FINLEY_DLL_API
545 virtual void setTags(const int functionSpaceType, const int newTag, const escript::Data& mask) const;
546
547 protected:
548
549 private:
550
551 //
552 // pointer to the externally created finley mesh
553 boost::shared_ptr<Finley_Mesh> m_finleyMesh;
554
555 static FunctionSpaceNamesMapType m_functionSpaceTypeNames;
556
557 };
558
559 } // end of namespace
560
561 #endif

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26