/[escript]/trunk/ripley/src/RipleyDomain.h
ViewVC logotype

Contents of /trunk/ripley/src/RipleyDomain.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3793 - (show annotations)
Wed Feb 1 07:39:43 2012 UTC (7 years, 8 months ago) by gross
File MIME type: text/plain
File size: 24510 byte(s)
new implementation of FCT solver with some modifications to the python interface
1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2012 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 #ifndef __RIPLEY_DOMAIN_H__
15 #define __RIPLEY_DOMAIN_H__
16
17 #include <ripley/Ripley.h>
18 #include <ripley/RipleyException.h>
19 #include <escript/AbstractContinuousDomain.h>
20 #include <escript/Data.h>
21 #include <escript/FunctionSpace.h>
22
23 struct Paso_Pattern;
24 struct Paso_SystemMatrixPattern;
25 struct Paso_SystemMatrix;
26
27 namespace ripley {
28
29 /**
30 \brief
31 RipleyDomain extends the AbstractContinuousDomain interface
32 for the Ripley library and is the base class for Rectangle and Brick.
33 */
34
35 class RipleyDomain : public escript::AbstractContinuousDomain
36 {
37 public:
38
39 /**
40 \brief recovers mesh from a file created with the dump() method
41 \param filename the name of the file
42 */
43 RIPLEY_DLL_API
44 static escript::Domain_ptr loadMesh(const std::string& filename);
45
46 /**
47 \brief reads a mesh from a file created with the write() method
48 \param filename the name of the file
49 */
50 RIPLEY_DLL_API
51 static escript::Domain_ptr readMesh(const std::string& filename);
52
53 /**
54 \brief
55 Constructor with number of dimensions. Allocates MPI info structure.
56 */
57 RIPLEY_DLL_API
58 RipleyDomain(dim_t dim);
59
60 /**
61 \brief
62 Destructor for RipleyDomain
63 */
64 RIPLEY_DLL_API
65 ~RipleyDomain();
66
67 /**
68 \brief
69 returns a description for this domain
70 */
71 RIPLEY_DLL_API
72 virtual std::string getDescription() const;
73
74 /**
75 \brief
76 returns the number of processors used for this domain
77 */
78 RIPLEY_DLL_API
79 virtual int getMPISize() const { return m_mpiInfo->size; }
80
81 /**
82 \brief
83 returns the MPI rank of this processor
84 */
85 RIPLEY_DLL_API
86 virtual int getMPIRank() const { return m_mpiInfo->rank; }
87
88 /**
89 \brief
90 if compiled for MPI then executes an MPI_Barrier, else does nothing
91 */
92 RIPLEY_DLL_API
93 virtual void MPIBarrier() const {
94 #ifdef ESYS_MPI
95 MPI_Barrier(m_mpiInfo->comm);
96 #endif
97 }
98
99 /**
100 \brief
101 returns true if on MPI processor 0, else false
102 */
103 RIPLEY_DLL_API
104 virtual bool onMasterProcessor() const { return getMPIRank()==0; }
105
106 /**
107 \brief
108 returns the MPI communicator
109 */
110 RIPLEY_DLL_API
111 #ifdef ESYS_MPI
112 MPI_Comm
113 #else
114 unsigned int
115 #endif
116 getMPIComm() const {
117 #ifdef ESYS_MPI
118 return m_mpiInfo->comm;
119 #else
120 return 0;
121 #endif
122 }
123
124 /**
125 \brief
126 returns true if the argument is a valid function space type for this
127 domain
128 */
129 RIPLEY_DLL_API
130 virtual bool isValidFunctionSpaceType(int fsType) const;
131
132 /**
133 \brief
134 returns a description for the given function space type code
135 */
136 RIPLEY_DLL_API
137 virtual std::string functionSpaceTypeAsString(int fsType) const;
138
139 /**
140 \brief
141 returns the number of spatial dimensions of the domain
142 */
143 RIPLEY_DLL_API
144 virtual int getDim() const { return m_numDim; }
145
146 /**
147 \brief equality operator
148 */
149 RIPLEY_DLL_API
150 virtual bool operator==(const escript::AbstractDomain& other) const;
151
152 /**
153 \brief inequality operator
154 */
155 RIPLEY_DLL_API
156 virtual bool operator!=(const escript::AbstractDomain& other) const {
157 return !(operator==(other));
158 }
159
160 /**
161 \brief
162 writes the current mesh to a file with the given name
163 \param filename The name of the file to write to
164 */
165 RIPLEY_DLL_API
166 void write(const std::string& filename) const;
167
168 /**
169 \brief
170 dumps the mesh to a file with the given name
171 \param filename The name of the output file
172 */
173 RIPLEY_DLL_API
174 void dump(const std::string& filename) const;
175
176 /**
177 \brief
178 returns the number of data points per sample, and the number of samples
179 as a pair.
180 \param fsType The function space type
181 */
182 RIPLEY_DLL_API
183 virtual std::pair<int,int> getDataShape(int fsType) const;
184
185 /**
186 \brief
187 returns the tag key for the given sample number
188 \param fsType The function space type
189 \param sampleNo The sample number
190 */
191 RIPLEY_DLL_API
192 int getTagFromSampleNo(int fsType, int sampleNo) const;
193
194 /**
195 \brief
196 sets a map from a clear tag name to a tag key
197 \param name tag name
198 \param tag tag key
199 */
200 RIPLEY_DLL_API
201 virtual void setTagMap(const std::string& name, int tag) {
202 m_tagMap[name] = tag;
203 }
204
205 /**
206 \brief
207 returns the tag key for tag name
208 \param name tag name
209 */
210 RIPLEY_DLL_API
211 virtual int getTag(const std::string& name) const {
212 if (m_tagMap.find(name) != m_tagMap.end()) {
213 return m_tagMap.find(name)->second;
214 } else {
215 throw RipleyException("getTag: invalid tag name");
216 }
217 }
218
219 /**
220 \brief
221 returns true if name is a defined tag name
222 \param name tag name to be checked
223 */
224 RIPLEY_DLL_API
225 virtual bool isValidTagName(const std::string& name) const {
226 return (m_tagMap.find(name)!=m_tagMap.end());
227 }
228
229 /**
230 \brief
231 returns all tag names in a single string separated by commas
232 */
233 RIPLEY_DLL_API
234 virtual std::string showTagNames() const;
235
236 /**
237 \brief
238 returns the array of reference numbers for a function space type
239 \param fsType The function space type
240 */
241 RIPLEY_DLL_API
242 const int* borrowSampleReferenceIDs(int fsType) const;
243
244 /**
245 \brief
246 assigns new location to the domain.
247 \note This is not supported in Ripley
248 */
249 RIPLEY_DLL_API
250 virtual void setNewX(const escript::Data& arg);
251
252 /**
253 \brief
254 interpolates data given on source onto target where source and target
255 have to be given on the same domain
256 */
257 RIPLEY_DLL_API
258 virtual void interpolateOnDomain(escript::Data& target, const escript::Data& source) const;
259
260 /**
261 \brief
262 returns true if data on fsType_source can be interpolated onto
263 fsType_target, false otherwise
264 */
265 RIPLEY_DLL_API
266 virtual bool probeInterpolationOnDomain(int fsType_source, int fsType_target) const;
267
268 /**
269 \brief
270 given a vector of FunctionSpace type codes, passes back a code which all
271 can be interpolated to
272 \return true if result is valid, false if not
273 */
274 RIPLEY_DLL_API
275 bool
276 commonFunctionSpace(const std::vector<int>& fs, int& resultcode) const;
277
278 /**
279 \brief
280 interpolates data given on source onto target where source and target
281 are given on different domains
282 */
283 RIPLEY_DLL_API
284 virtual void interpolateACross(escript::Data& target, const escript::Data& source) const;
285
286 /**
287 \brief
288 determines whether interpolation from source to target is possible
289 */
290 RIPLEY_DLL_API
291 virtual bool probeInterpolationACross(int, const escript::AbstractDomain&, int) const;
292
293 /**
294 \brief
295 returns locations in the FEM nodes
296 */
297 RIPLEY_DLL_API
298 virtual escript::Data getX() const;
299
300 /**
301 \brief
302 returns boundary normals at the quadrature point on the face elements
303 */
304 RIPLEY_DLL_API
305 virtual escript::Data getNormal() const;
306
307 /**
308 \brief returns the element size
309 */
310 RIPLEY_DLL_API
311 virtual escript::Data getSize() const;
312
313 /**
314 \brief
315 copies the location of data points into arg. The domain of arg has to
316 match this domain.
317 */
318 RIPLEY_DLL_API
319 virtual void setToX(escript::Data& arg) const;
320
321 /**
322 \brief
323 copies the surface normals at data points into out. The actual function
324 space to be considered is defined by out. out has to be defined on this
325 domain.
326 */
327 RIPLEY_DLL_API
328 virtual void setToNormal(escript::Data& out) const;
329
330 /**
331 \brief
332 copies the size of samples into out. The actual function space to be
333 considered is defined by out. out has to be defined on this domain.
334 */
335 RIPLEY_DLL_API
336 virtual void setToSize(escript::Data& out) const;
337
338 /**
339 \brief
340 copies the gradient of 'in' into 'out'. The actual function space to be
341 considered for the gradient is defined by 'in'. Both arguments have to
342 be defined on this domain.
343 */
344 RIPLEY_DLL_API
345 virtual void setToGradient(escript::Data& out, const escript::Data& in) const;
346
347 /**
348 \brief
349 returns true if this rank owns the sample id on given function space
350 */
351 RIPLEY_DLL_API
352 virtual bool ownSample(int fsType, index_t id) const;
353
354 /**
355 \brief
356 returns the number of data points summed across all MPI processes
357 */
358 RIPLEY_DLL_API
359 virtual int getNumDataPointsGlobal() const;
360
361 /**
362 \brief
363 assigns new tag newTag to all samples of given function space with a
364 positive value of mask for any of its sample points
365 */
366 RIPLEY_DLL_API
367 virtual void setTags(const int fsType, const int newTag, const escript::Data& mask) const;
368
369 /**
370 \brief
371 returns true if data on this domain and given function space type has
372 to be considered as cell centered data
373 */
374 RIPLEY_DLL_API
375 virtual bool isCellOriented(int fsType) const;
376
377 /**
378 \brief
379 returns a status indicator of the domain. The status identifier should
380 be unique over the lifetime of the object but may be updated if changes
381 to the domain happen, e.g. modifications to its geometry.
382 */
383 RIPLEY_DLL_API
384 virtual StatusType getStatus() const { return m_status; }
385
386 /**
387 \brief
388 returns the number of tags in use for a function space type
389 */
390 RIPLEY_DLL_API
391 virtual int getNumberOfTagsInUse(int fsType) const;
392
393 /**
394 \brief
395 returns a pointer to the list of tags in use for a function space type
396 */
397 RIPLEY_DLL_API
398 virtual const int* borrowListOfTagsInUse(int fsType) const;
399
400 /**
401 \brief
402 checks if this domain allows tags for the specified function space type
403 */
404 RIPLEY_DLL_API
405 virtual bool canTag(int fsType) const;
406
407 /**
408 \brief
409 returns the approximation order used for a function space
410 */
411 RIPLEY_DLL_API
412 virtual int getApproximationOrder(const int fsType) const { return 1; }
413
414 /**
415 \brief
416 returns true if this domain supports contact elements, false otherwise
417 */
418 RIPLEY_DLL_API
419 virtual bool supportsContactElements() const { return false; }
420
421 /**
422 \brief
423 returns a continuous FunctionSpace code
424 */
425 RIPLEY_DLL_API
426 virtual int getContinuousFunctionCode() const { return Nodes; }
427
428 /**
429 \brief
430 returns a continuous on reduced order nodes FunctionSpace code
431 */
432 RIPLEY_DLL_API
433 virtual int getReducedContinuousFunctionCode() const { return ReducedNodes; }
434
435 /**
436 \brief
437 returns a function FunctionSpace code
438 */
439 RIPLEY_DLL_API
440 virtual int getFunctionCode() const { return Elements; }
441
442 /**
443 \brief
444 returns a function with reduced integration order FunctionSpace code
445 */
446 RIPLEY_DLL_API
447 virtual int getReducedFunctionCode() const { return ReducedElements; }
448
449 /**
450 \brief
451 returns a function on boundary FunctionSpace code
452 */
453 RIPLEY_DLL_API
454 virtual int getFunctionOnBoundaryCode() const { return FaceElements; }
455
456 /**
457 \brief
458 returns a function on boundary with reduced integration order
459 FunctionSpace code
460 */
461 RIPLEY_DLL_API
462 virtual int getReducedFunctionOnBoundaryCode() const { return ReducedFaceElements; }
463
464 /**
465 \brief
466 return a FunctionOnContactZero code
467 */
468 RIPLEY_DLL_API
469 virtual int getFunctionOnContactZeroCode() const {
470 throw RipleyException("Ripley does not support contact elements");
471 }
472
473 /**
474 \brief
475 returns a FunctionOnContactZero code with reduced integration order
476 */
477 RIPLEY_DLL_API
478 virtual int getReducedFunctionOnContactZeroCode() const {
479 throw RipleyException("Ripley does not support contact elements");
480 }
481
482 /**
483 \brief
484 returns a FunctionOnContactOne code
485 */
486 RIPLEY_DLL_API
487 virtual int getFunctionOnContactOneCode() const {
488 throw RipleyException("Ripley does not support contact elements");
489 }
490
491 /**
492 \brief
493 returns a FunctionOnContactOne code with reduced integration order
494 */
495 RIPLEY_DLL_API
496 virtual int getReducedFunctionOnContactOneCode() const {
497 throw RipleyException("Ripley does not support contact elements");
498 }
499
500 /**
501 \brief
502 returns a Solution FunctionSpace code
503 */
504 RIPLEY_DLL_API
505 virtual int getSolutionCode() const { return DegreesOfFreedom; }
506
507 /**
508 \brief
509 returns a ReducedSolution FunctionSpace code
510 */
511 RIPLEY_DLL_API
512 virtual int getReducedSolutionCode() const { return ReducedDegreesOfFreedom; }
513
514 /**
515 \brief
516 returns a DiracDeltaFunctions FunctionSpace code
517 */
518 RIPLEY_DLL_API
519 virtual int getDiracDeltaFunctionsCode() const { return Points; }
520
521 /**
522 \brief
523 returns the identifier of the matrix type to be used for the global
524 stiffness matrix when a particular solver, package, preconditioner,
525 and symmetric matrix is used
526 \param solver
527 \param preconditioner
528 \param package
529 \param symmetry
530 */
531 RIPLEY_DLL_API
532 virtual int getSystemMatrixTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
533
534 /**
535 \brief
536 returns the identifier of the transport problem type to be used when a
537 particular solver, preconditioner, package and symmetric matrix is used
538 \param solver
539 \param preconditioner
540 \param package
541 \param symmetry
542 */
543 RIPLEY_DLL_API
544 virtual int getTransportTypeId(const int solver, const int preconditioner, const int package, const bool symmetry) const;
545
546 /**
547 \brief
548 copies the integrals of the function defined by arg into integrals.
549 arg has to be defined on this domain.
550 */
551 RIPLEY_DLL_API
552 virtual void setToIntegrals(std::vector<double>& integrals, const escript::Data& arg) const;
553
554 /**
555 \brief
556 adds a PDE onto the stiffness matrix mat and rhs
557 */
558 RIPLEY_DLL_API
559 virtual void addPDEToSystem(escript::AbstractSystemMatrix& mat,
560 escript::Data& rhs, const escript::Data& A, const escript::Data& B,
561 const escript::Data& C, const escript::Data& D,
562 const escript::Data& X, const escript::Data& Y,
563 const escript::Data& d, const escript::Data& y,
564 const escript::Data& d_contact, const escript::Data& y_contact,
565 const escript::Data& d_dirac, const escript::Data& y_dirac) const;
566
567 /**
568 \brief
569 adds a PDE onto rhs
570 */
571 RIPLEY_DLL_API
572 virtual void addPDEToRHS(escript::Data& rhs, const escript::Data& X,
573 const escript::Data& Y, const escript::Data& y,
574 const escript::Data& y_contact, const escript::Data& y_dirac) const;
575
576 /**
577 \brief
578 adds a PDE onto a transport problem
579 */
580 RIPLEY_DLL_API
581 virtual void addPDEToTransportProblem(escript::AbstractTransportProblem& tp,
582 escript::Data& source, const escript::Data& M,
583 const escript::Data& A, const escript::Data& B,
584 const escript::Data& C, const escript::Data& D,
585 const escript::Data& X, const escript::Data& Y,
586 const escript::Data& d, const escript::Data& y,
587 const escript::Data& d_contact, const escript::Data& y_contact,
588 const escript::Data& d_dirac, const escript::Data& y_dirac) const;
589
590
591 /**
592 \brief
593 creates a stiffness matrix and initializes it with zeros
594 */
595 RIPLEY_DLL_API
596 virtual escript::ASM_ptr newSystemMatrix(const int row_blocksize,
597 const escript::FunctionSpace& row_functionspace,
598 const int column_blocksize,
599 const escript::FunctionSpace& column_functionspace, const int type) const;
600
601 /**
602 \brief
603 creates a transport problem
604 */
605 RIPLEY_DLL_API
606 virtual escript::ATP_ptr newTransportProblem(
607 const int blocksize, const escript::FunctionSpace& functionspace,
608 const int type) const;
609
610 /**
611 \brief
612 writes information about the mesh to standard output
613 \param full whether to print additional data
614 */
615 RIPLEY_DLL_API
616 virtual void Print_Mesh_Info(const bool full=false) const;
617
618 /**
619 \brief
620 returns the number of nodes per MPI rank in each dimension
621 */
622 RIPLEY_DLL_API
623 virtual IndexVector getNumNodesPerDim() const;
624
625 /**
626 \brief
627 returns the number of elements per MPI rank in each dimension
628 */
629 RIPLEY_DLL_API
630 virtual IndexVector getNumElementsPerDim() const;
631
632 /**
633 \brief
634 returns the number of face elements in the order
635 (left,right,bottom,top,[front,back]) on current MPI rank
636 */
637 RIPLEY_DLL_API
638 virtual IndexVector getNumFacesPerBoundary() const;
639
640 /**
641 \brief
642 returns the node distribution vector
643 */
644 RIPLEY_DLL_API
645 virtual IndexVector getNodeDistribution() const;
646
647 /**
648 \brief
649 returns the number of spatial subdivisions in each dimension
650 */
651 RIPLEY_DLL_API
652 virtual IndexVector getNumSubdivisionsPerDim() const;
653
654 /**
655 \brief
656 returns the first coordinate value and the node spacing along given
657 dimension as a pair
658 */
659 RIPLEY_DLL_API
660 virtual std::pair<double,double> getFirstCoordAndSpacing(dim_t dim) const;
661
662 protected:
663 /// returns the number of nodes per MPI rank
664 virtual dim_t getNumNodes() const;
665
666 /// returns the number of elements per MPI rank
667 virtual dim_t getNumElements() const;
668
669 /// returns the number of degrees of freedom per MPI rank
670 virtual dim_t getNumDOF() const;
671
672 /// returns the number of face elements on current MPI rank
673 virtual dim_t getNumFaceElements() const;
674
675 /// inserts the nodes that share an element with 'node' into 'index' and
676 /// returns the number of these neighbours
677 virtual dim_t insertNeighbourNodes(IndexVector& index, index_t node) const;
678
679 /// populates the data object 'arg' with the node coordinates
680 virtual void assembleCoordinates(escript::Data& arg) const;
681
682 /// computes the gradient of 'in' and puts the result in 'out'
683 virtual void assembleGradient(escript::Data& out, escript::Data& in) const;
684
685 /// copies the integrals of the function defined by 'arg' into 'integrals'
686 virtual void assembleIntegrate(std::vector<double>& integrals, escript::Data& arg) const;
687
688 /// assembles a single PDE into the system matrix 'mat' and the right hand
689 /// side 'rhs'
690 virtual void assemblePDESingle(Paso_SystemMatrix* mat, escript::Data& rhs,
691 const escript::Data& A, const escript::Data& B,
692 const escript::Data& C, const escript::Data& D,
693 const escript::Data& X, const escript::Data& Y) const;
694
695 /// assembles boundary conditions of a single PDE into the system matrix
696 /// 'mat' and the right hand side 'rhs'
697 virtual void assemblePDEBoundarySingle(Paso_SystemMatrix* mat,
698 escript::Data& rhs, const escript::Data& d,
699 const escript::Data& y) const;
700
701 /// assembles a single PDE with reduced order into the system matrix 'mat'
702 /// and the right hand side 'rhs'
703 virtual void assemblePDESingleReduced(Paso_SystemMatrix* mat,
704 escript::Data& rhs, const escript::Data& A, const escript::Data& B,
705 const escript::Data& C, const escript::Data& D,
706 const escript::Data& X, const escript::Data& Y) const;
707
708 /// assembles boundary conditions of a single PDE with reduced order into
709 /// the system matrix 'mat' and the right hand side 'rhs'
710 virtual void assemblePDEBoundarySingleReduced(Paso_SystemMatrix* mat,
711 escript::Data& rhs, const escript::Data& d,
712 const escript::Data& y) const;
713
714 /// assembles a system of PDEs into the system matrix 'mat' and the right
715 /// hand side 'rhs'
716 virtual void assemblePDESystem(Paso_SystemMatrix* mat, escript::Data& rhs,
717 const escript::Data& A, const escript::Data& B,
718 const escript::Data& C, const escript::Data& D,
719 const escript::Data& X, const escript::Data& Y) const;
720
721 /// assembles boundary conditions of a system of PDEs into the system
722 /// matrix 'mat' and the right hand side 'rhs'
723 virtual void assemblePDEBoundarySystem(Paso_SystemMatrix* mat,
724 escript::Data& rhs, const escript::Data& d,
725 const escript::Data& y) const;
726
727 /// assembles a system of PDEs with reduced order into the system matrix
728 /// 'mat' and the right hand side 'rhs'
729 virtual void assemblePDESystemReduced(Paso_SystemMatrix* mat,
730 escript::Data& rhs, const escript::Data& A, const escript::Data& B,
731 const escript::Data& C, const escript::Data& D,
732 const escript::Data& X, const escript::Data& Y) const;
733
734 /// assembles boundary conditions of a system of PDEs with reduced order
735 /// into the system matrix 'mat' and the right hand side 'rhs'
736 virtual void assemblePDEBoundarySystemReduced(Paso_SystemMatrix* mat,
737 escript::Data& rhs, const escript::Data& d,
738 const escript::Data& y) const;
739
740 /// returns the Paso system matrix pattern
741 virtual Paso_SystemMatrixPattern* getPattern(bool reducedRowOrder,
742 bool reducedColOrder) const;
743
744 /// interpolates data on nodes in 'in' onto (reduced) elements in 'out'
745 virtual void interpolateNodesOnElements(escript::Data& out,
746 escript::Data& in, bool reduced) const;
747
748 /// interpolates data on nodes in 'in' onto (reduced) face elements in 'out'
749 virtual void interpolateNodesOnFaces(escript::Data& out, escript::Data& in,
750 bool reduced) const;
751
752 /// converts data on nodes in 'in' to degrees of freedom in 'out'
753 virtual void nodesToDOF(escript::Data& out, escript::Data& in) const;
754
755 /// converts data on degrees of freedom in 'in' to nodes in 'out'
756 virtual void dofToNodes(escript::Data& out, escript::Data& in) const;
757
758 /// copies data in 'in' to 'out' (both must be on same function space)
759 void copyData(escript::Data& out, escript::Data& in) const;
760
761 /// averages data in 'in' to 'out' (from non-reduced to reduced fs)
762 void averageData(escript::Data& out, escript::Data& in) const;
763
764 // this is const because setTags is const
765 void updateTagsInUse(int fsType) const;
766
767 /// allocates and returns a Paso pattern structure
768 Paso_Pattern* createPasoPattern(const IndexVector& ptr,
769 const IndexVector& index, const dim_t M, const dim_t N) const;
770
771 /// creates the pattern for the main block of the system matrix
772 Paso_Pattern* createMainPattern() const;
773
774 /// creates the pattern for the column and row couple blocks of the system
775 /// matrix. colIndices[i] contains all IDs of DOFs that are connected with
776 /// DOF i but remote and 'N' is the total number of remote components
777 void createCouplePatterns(const std::vector<IndexVector>& colIndices,
778 const dim_t N, Paso_Pattern** colPattern,
779 Paso_Pattern** rowPattern) const;
780
781 void addToSystemMatrix(Paso_SystemMatrix* in, const IndexVector& nodes_Eq,
782 dim_t num_Eq, const IndexVector& nodes_Sol, dim_t num_Sol,
783 const std::vector<double>& array) const;
784
785 dim_t m_numDim;
786 StatusType m_status;
787 Esys_MPIInfo *m_mpiInfo;
788 TagMap m_tagMap;
789 mutable IndexVector m_nodeTags, m_nodeTagsInUse;
790 mutable IndexVector m_elementTags, m_elementTagsInUse;
791 mutable IndexVector m_faceTags, m_faceTagsInUse;
792 };
793
794 } // end of namespace ripley
795
796 #endif // __RIPLEY_DOMAIN_H__
797

  ViewVC Help
Powered by ViewVC 1.1.26