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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4657 - (show annotations)
Thu Feb 6 06:12:20 2014 UTC (5 years, 7 months ago) by jfenwick
File MIME type: text/plain
File size: 25778 byte(s)
I changed some files.
Updated copyright notices, added GeoComp.



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

  ViewVC Help
Powered by ViewVC 1.1.26