/[escript]/branches/diaplayground/ripley/src/RipleyDomain.h
ViewVC logotype

Contents of /branches/diaplayground/ripley/src/RipleyDomain.h

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5136 - (show annotations)
Tue Sep 9 07:13:55 2014 UTC (4 years, 5 months ago) by caltinay
File MIME type: text/plain
File size: 25338 byte(s)
ripley now supports paso solvers again and returns an appropriate matrix type
id. Changed the getSystemMatrixTypeId() method to take a full SolverBuddy
instance and made some other simplifications.

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

  ViewVC Help
Powered by ViewVC 1.1.26