/[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 2521 - (show annotations)
Tue Jul 7 00:08:58 2009 UTC (9 years, 9 months ago) by jfenwick
File MIME type: text/plain
File size: 17865 byte(s)
Modified Lazy so that resolving a single sample uses the node cache method.
Fixed some doxygen problems.
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2008 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14
15 #if !defined finley_MeshAdapter_20040526_H
16 #define finley_MeshAdapter_20040526_H
17 #include "system_dep.h"
18
19 extern "C" {
20 #include "../Mesh.h"
21 #include "../Finley.h"
22 #include "../Assemble.h"
23 #include "paso/SystemMatrix.h"
24 #include "paso/SolverFCT.h"
25 #include "paso/Paso_MPI.h"
26 }
27
28 #include "FinleyError.h"
29 #include "FinleyAdapterException.h"
30
31 #include "SystemMatrixAdapter.h"
32 #include "TransportProblemAdapter.h"
33 #include "escript/AbstractContinuousDomain.h"
34 #include "escript/FunctionSpace.h"
35 #include "escript/FunctionSpaceFactory.h"
36
37 #include <boost/shared_ptr.hpp>
38 #include <boost/python/dict.hpp>
39 #include <boost/python/extract.hpp>
40
41 #include <map>
42 #include <vector>
43 #include <string>
44 #include <sstream>
45
46 //
47 // forward declarations
48 class Data;
49
50 //using namespace escript;
51
52 namespace finley {
53
54 struct null_deleter
55 {
56 void operator()(void const *ptr) const
57 {
58 }
59 };
60
61
62 /**
63 \brief
64 MeshAdapter implements the AbstractContinuousDomain
65 interface for the Finley library.
66
67 Description:
68 MeshAdapter implements the AbstractContinuousDomain
69 interface for the Finley library.
70 */
71
72 class MeshAdapter : public escript::AbstractContinuousDomain {
73
74 public:
75
76 //
77 // Codes for function space types supported
78 static const int DegreesOfFreedom;
79 static const int ReducedDegreesOfFreedom;
80 static const int Nodes;
81 static const int ReducedNodes;
82 static const int Elements;
83 static const int ReducedElements;
84 static const int FaceElements;
85 static const int ReducedFaceElements;
86 static const int Points;
87 static const int ContactElementsZero;
88 static const int ReducedContactElementsZero;
89 static const int ContactElementsOne;
90 static const int ReducedContactElementsOne;
91
92 /**
93 \brief
94 Constructor for MeshAdapter
95
96 Description:
97 Constructor for MeshAdapter. The pointer passed to MeshAdapter
98 is deleted using a call to Finley_Mesh_free in the
99 MeshAdapter destructor.
100
101 Throws:
102 May throw an exception derived from EsysException
103
104 \param finleyMesh Input - A pointer to the externally constructed
105 finley mesh.The pointer passed to MeshAdapter
106 is deleted using a call to
107 Finley_Mesh_free in the MeshAdapter
108 destructor.
109 */
110 FINLEY_DLL_API
111 MeshAdapter(Finley_Mesh* finleyMesh=0);
112
113 /**
114 \brief
115 Copy constructor.
116 */
117 FINLEY_DLL_API
118 MeshAdapter(const MeshAdapter& in);
119
120 /**
121 \brief
122 Destructor for MeshAdapter. As specified in the constructor
123 this calls Finley_Mesh_free for the pointer given to the
124 constructor.
125 */
126 FINLEY_DLL_API
127 ~MeshAdapter();
128
129 /**
130 \brief
131 return the number of processors used for this domain
132 */
133 FINLEY_DLL_API
134 virtual int getMPISize() const;
135 /**
136 \brief
137 return the number MPI rank of this processor
138 */
139
140 FINLEY_DLL_API
141 virtual int getMPIRank() const;
142
143 /**
144 \brief
145 If compiled for MPI then execute an MPI_Barrier, else do nothing
146 */
147
148 FINLEY_DLL_API
149 virtual void MPIBarrier() const;
150
151 /**
152 \brief
153 Return true if on MPI processor 0, else false
154 */
155
156 FINLEY_DLL_API
157 virtual bool onMasterProcessor() const;
158
159 /**
160 \brief
161 return this as an AbstractContinuousDomain.
162 */
163 inline const AbstractContinuousDomain& asAbstractContinuousDomain() const
164 {
165 return *(static_cast<const AbstractContinuousDomain*>(this));
166 }
167
168 /**
169 \brief
170 return this as an AbstractDomain.
171 */
172 inline const AbstractDomain& asAbstractDomain() const
173 {
174 return *(static_cast<const AbstractDomain*>(this));
175 }
176
177 /**
178 \brief
179 Write the current mesh to a file with the given name.
180 \param fileName Input - The name of the file to write to.
181 */
182 FINLEY_DLL_API
183 void write(const std::string& fileName) const;
184
185 /**
186 \brief
187 \param full
188 */
189 FINLEY_DLL_API
190 void Print_Mesh_Info(const bool full=false) const;
191
192 /**
193 \brief
194 dumps the mesh to a file with the given name.
195 \param fileName Input - The name of the file
196 */
197 FINLEY_DLL_API
198 void dump(const std::string& fileName) const;
199
200 /**
201 \brief
202 return the pointer to the underlying finley mesh structure
203 */
204 FINLEY_DLL_API
205 Finley_Mesh* getFinley_Mesh() const;
206
207 /**
208 \brief
209 Return the tag key for the given sample number.
210 \param functionSpaceType Input - The function space type.
211 \param sampleNo Input - The sample number.
212 */
213 FINLEY_DLL_API
214 int getTagFromSampleNo(int functionSpaceType, int sampleNo) const;
215
216 /**
217 \brief
218 Return the reference number of the given sample number.
219 \param functionSpaceType Input - The function space type.
220 */
221 FINLEY_DLL_API
222 const int* borrowSampleReferenceIDs(int functionSpaceType) const;
223
224 /**
225 \brief
226 Returns true if the given integer is a valid function space type
227 for this domain.
228 */
229 FINLEY_DLL_API
230 virtual bool isValidFunctionSpaceType(int functionSpaceType) const;
231
232 /**
233 \brief
234 Return a description for this domain
235 */
236 FINLEY_DLL_API
237 virtual std::string getDescription() const;
238
239 /**
240 \brief
241 Return a description for the given function space type code
242 */
243 FINLEY_DLL_API
244 virtual std::string functionSpaceTypeAsString(int functionSpaceType) const;
245
246 /**
247 \brief
248 Build the table of function space type names
249 */
250 FINLEY_DLL_API
251 void setFunctionSpaceTypeNames();
252
253 /**
254 \brief
255 Return a continuous FunctionSpace code
256 */
257 FINLEY_DLL_API
258 virtual int getContinuousFunctionCode() const;
259
260 /**
261 \brief
262 Return a continuous on reduced order nodes FunctionSpace code
263 */
264 FINLEY_DLL_API
265 virtual int getReducedContinuousFunctionCode() const;
266
267 /**
268 \brief
269 Return a function FunctionSpace code
270 */
271 FINLEY_DLL_API
272 virtual int getFunctionCode() const;
273
274 /**
275 \brief
276 Return a function with reduced integration order FunctionSpace code
277 */
278 FINLEY_DLL_API
279 virtual int getReducedFunctionCode() const;
280
281 /**
282 \brief
283 Return a function on boundary FunctionSpace code
284 */
285 FINLEY_DLL_API
286 virtual int getFunctionOnBoundaryCode() const;
287
288 /**
289 \brief
290 Return a function on boundary with reduced integration order FunctionSpace code
291 */
292 FINLEY_DLL_API
293 virtual int getReducedFunctionOnBoundaryCode() const;
294
295 /**
296 \brief
297 Return a FunctionOnContactZero code
298 */
299 FINLEY_DLL_API
300 virtual int getFunctionOnContactZeroCode() const;
301
302 /**
303 \brief
304 Return a FunctionOnContactZero code with reduced integration order
305 */
306 FINLEY_DLL_API
307 virtual int getReducedFunctionOnContactZeroCode() const;
308
309 /**
310 \brief
311 Return a FunctionOnContactOne code
312 */
313 FINLEY_DLL_API
314 virtual int getFunctionOnContactOneCode() const;
315
316 /**
317 \brief
318 Return a FunctionOnContactOne code with reduced integration order
319 */
320 FINLEY_DLL_API
321 virtual int getReducedFunctionOnContactOneCode() const;
322
323 /**
324 \brief
325 Return a Solution code
326 */
327 FINLEY_DLL_API
328 virtual int getSolutionCode() const;
329
330 /**
331 \brief
332 Return a ReducedSolution code
333 */
334 FINLEY_DLL_API
335 virtual int getReducedSolutionCode() const;
336
337 /**
338 \brief
339 Return a DiracDeltaFunction code
340 */
341 FINLEY_DLL_API
342 virtual int getDiracDeltaFunctionCode() const;
343
344 /**
345 5B
346 \brief
347 */
348 typedef std::map<int, std::string> FunctionSpaceNamesMapType;
349
350 /**
351 \brief
352 */
353 FINLEY_DLL_API
354 virtual int getDim() const;
355
356 /**
357 \brief
358 Return the number of data points summed across all MPI processes
359 */
360 FINLEY_DLL_API
361 virtual int getNumDataPointsGlobal() const;
362
363 /**
364 \brief
365 Return the number of data points per sample, and the number of samples as a pair.
366 \param functionSpaceCode Input -
367 */
368 FINLEY_DLL_API
369 virtual std::pair<int,int> getDataShape(int functionSpaceCode) const;
370
371 /**
372 \brief
373 copies the location of data points into arg. The domain of arg has to match this.
374 has to be implemented by the actual Domain adapter.
375 */
376 FINLEY_DLL_API
377 virtual void setToX(escript::Data& arg) const;
378
379 /**
380 \brief
381 sets a map from a clear tag name to a tag key
382 \param name Input - tag name.
383 \param tag Input - tag key.
384 */
385 FINLEY_DLL_API
386 virtual void setTagMap(const std::string& name, int tag);
387
388 /**
389 \brief
390 Return the tag key for tag name.
391 \param name Input - tag name
392 */
393 FINLEY_DLL_API
394 virtual int getTag(const std::string& name) const;
395
396 /**
397 \brief
398 Returns true if name is a defined tage name.
399 \param name Input - tag name to be checked.
400 */
401 FINLEY_DLL_API
402 virtual bool isValidTagName(const std::string& name) const;
403
404 /**
405 \brief
406 Returns all tag names in a single string sperated by commas
407 */
408 FINLEY_DLL_API
409 virtual std::string showTagNames() const;
410
411 /**
412 \brief
413 assigns new location to the domain
414 */
415 FINLEY_DLL_API
416 virtual void setNewX(const escript::Data& arg);
417
418 /**
419 \brief
420 interpolates data given on source onto target where source and target have to be given on the same domain.
421 */
422 FINLEY_DLL_API
423 virtual void interpolateOnDomain(escript::Data& target,const escript::Data& source) const;
424
425
426 FINLEY_DLL_API
427 virtual bool probeInterpolationOnDomain(int functionSpaceType_source,int functionSpaceType_target) const;
428
429 /**
430 \brief
431 interpolates data given on source onto target where source and target are given on different domains.
432 has to be implemented by the actual Domain adapter.
433 */
434 FINLEY_DLL_API
435 virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
436
437 /**
438 \brief determines whether interpolation from source to target is possible.
439 Must be implemented by the actual Domain adapter
440 */
441 FINLEY_DLL_API
442 virtual bool probeInterpolationACross(int functionSpaceType_source,const escript::AbstractDomain& targetDomain, int functionSpaceType_target) const;
443
444 /**
445 \brief
446 copies the surface normals at data points into out. The actual function space to be considered
447 is defined by out. out has to be defined on this.
448 */
449 FINLEY_DLL_API
450 virtual void setToNormal(escript::Data& out) const;
451
452 /**
453 \brief
454 copies the size of samples into out. The actual function space to be considered
455 is defined by out. out has to be defined on this.
456 */
457 FINLEY_DLL_API
458 virtual void setToSize(escript::Data& out) const;
459
460 /**
461 \brief
462 copies the gradient of arg into grad. The actual function space to be considered
463 for the gradient is defined by grad. arg and grad have to be defined on this.
464 */
465 FINLEY_DLL_API
466 virtual void setToGradient(escript::Data& grad,const escript::Data& arg) const;
467
468 /**
469 \brief
470 copies the integrals of the function defined by arg into integrals.
471 arg has to be defined on this.
472 */
473 FINLEY_DLL_API
474 virtual void setToIntegrals(std::vector<double>& integrals,const escript::Data& arg) const;
475
476 /**
477 \brief
478 return the identifier of the matrix type to be used for the global stiffness matrix when a particular solver, package, perconditioner,
479 and symmetric matrix is used.
480 \param solver
481 \param preconditioner
482 \param package
483 \param symmetry
484 */
485 FINLEY_DLL_API
486 virtual int getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
487
488 /**
489 \brief
490 return the identifier of the transport problem type to be used when a particular solver, perconditioner, package
491 and symmetric matrix is used.
492 \param solver
493 \param preconditioner
494 \param package
495 \param symmetry
496 */
497 FINLEY_DLL_API
498 virtual int getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
499
500 /**
501 \brief
502 returns true if data on this domain and a function space of type functionSpaceCode has to
503 considered as cell centered data.
504 */
505 FINLEY_DLL_API
506 virtual bool isCellOriented(int functionSpaceCode) const;
507
508 /**
509 \brief
510 Saves a dictonary of Data objects to an OpenDX input file. The keywords are used as identifier
511
512 This has to be implemented by the actual Domain adapter.
513 */
514 FINLEY_DLL_API
515 virtual void saveDX(const std::string& filename,const boost::python::dict& arg) const;
516
517
518 /**
519 \brief
520 Saves a dictonary of Data objects to an VTK XML input file. The keywords are used as identifier
521
522 This has to be implemented by the actual Domain adapter.
523 */
524 FINLEY_DLL_API
525 virtual void saveVTK(const std::string& filename,const boost::python::dict& arg, const std::string& metadata, const std::string& metadata_schema) const;
526
527 /**
528 \brief
529 returns the function space representation of the type functionSpaceCode on this domain
530 as a vtkObject.
531 */
532 // vtkObject createVtkObject(int functionSpaceCode) const;
533
534 /**
535 \brief
536 adds a PDE onto the stiffness matrix mat and a rhs
537 */
538 FINLEY_DLL_API
539 virtual void addPDEToSystem(
540 SystemMatrixAdapter& mat, escript::Data& rhs,
541 const escript::Data& A, const escript::Data& B, const escript::Data& C,
542 const escript::Data& D, const escript::Data& X, const escript::Data& Y,
543 const escript::Data& d, const escript::Data& y,
544 const escript::Data& d_contact, const escript::Data& y_contact) const;
545 /**
546 \brief
547 adds a PDE onto the lumped stiffness matrix matrix
548 */
549 FINLEY_DLL_API
550 virtual void addPDEToLumpedSystem(
551 escript::Data& mat,
552 const escript::Data& D,
553 const escript::Data& d) const;
554
555 /**
556 \brief
557 adds a PDE onto the stiffness matrix mat and a rhs
558 */
559 FINLEY_DLL_API
560 virtual void addPDEToRHS(escript::Data& rhs,
561 const escript::Data& X, const escript::Data& Y,
562 const escript::Data& y, const escript::Data& y_contact) const;
563 /**
564 \brief
565 adds a PDE onto a transport problem
566 */
567
568 FINLEY_DLL_API
569 virtual void addPDEToTransportProblem(
570 TransportProblemAdapter& tp, escript::Data& source,
571 const escript::Data& M,
572 const escript::Data& A, const escript::Data& B, const escript::Data& C,const escript::Data& D,
573 const escript::Data& X,const escript::Data& Y,
574 const escript::Data& d, const escript::Data& y,
575 const escript::Data& d_contact,const escript::Data& y_contact) const;
576
577
578 /**
579 \brief
580 creates a SystemMatrixAdapter stiffness matrix and initializes it with zeros:
581 */
582 FINLEY_DLL_API
583 SystemMatrixAdapter newSystemMatrix(
584 const int row_blocksize,
585 const escript::FunctionSpace& row_functionspace,
586 const int column_blocksize,
587 const escript::FunctionSpace& column_functionspace,
588 const int type) const;
589 /**
590 \brief
591 creates a TransportProblemAdapter
592
593 */
594
595 FINLEY_DLL_API
596 TransportProblemAdapter newTransportProblem(
597 const double theta,
598 const int blocksize,
599 const escript::FunctionSpace& functionspace,
600 const int type) const;
601
602 /**
603 \brief returns locations in the FEM nodes
604 */
605 FINLEY_DLL_API
606 virtual escript::Data getX() const;
607
608 /**
609 \brief return boundary normals at the quadrature point on the face elements
610 */
611 FINLEY_DLL_API
612 virtual escript::Data getNormal() const;
613
614 /**
615 \brief returns the element size
616 */
617 FINLEY_DLL_API
618 virtual escript::Data getSize() const;
619
620 /**
621 \brief comparison operators
622 */
623 FINLEY_DLL_API
624 virtual bool operator==(const escript::AbstractDomain& other) const;
625 FINLEY_DLL_API
626 virtual bool operator!=(const escript::AbstractDomain& other) const;
627
628 /**
629 \brief assigns new tag newTag to all samples of functionspace with a positive
630 value of mask for any its sample point.
631
632 */
633 FINLEY_DLL_API
634 virtual void setTags(const int functionSpaceType, const int newTag, const escript::Data& mask) const;
635
636 /**
637 \brief
638 return the number of tags in use and a pointer to an array with the number of tags in use
639 */
640 FINLEY_DLL_API
641 virtual int getNumberOfTagsInUse(int functionSpaceCode) const;
642
643 FINLEY_DLL_API
644 virtual const int* borrowListOfTagsInUse(int functionSpaceCode) const;
645
646
647 /**
648 \brief Checks if this domain allows tags for the specified functionSpaceCode.
649 */
650 FINLEY_DLL_API
651 virtual
652 bool canTag(int functionSpaceCode) const;
653
654
655 protected:
656
657 private:
658 void extractArgsFromDict(const boost::python::dict& arg, int& numData,
659 char**& names, escriptDataC*& data,
660 escriptDataC**& dataPtr) const;
661
662 //
663 // pointer to the externally created finley mesh
664 boost::shared_ptr<Finley_Mesh> m_finleyMesh;
665
666 static FunctionSpaceNamesMapType m_functionSpaceTypeNames;
667
668 };
669
670 } // end of namespace
671
672 #endif

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26