libMesh
petsc_macro.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 #ifndef LIBMESH_PETSC_MACRO_H
19 #define LIBMESH_PETSC_MACRO_H
20 
21 // Local includes
22 #include "libmesh/libmesh_config.h"
23 
24 #ifdef LIBMESH_HAVE_PETSC
25 
26 #include <petscversion.h>
27 
28 // A convenient macro for comparing PETSc versions. Returns 1 if the
29 // current PETSc version is < major.minor.subminor and zero otherwise.
30 #define PETSC_VERSION_LESS_THAN(major,minor,subminor) \
31  PETSC_VERSION_LT(major,minor,subminor)
32 
33 #define PETSC_VERSION_EQUALS(major,minor,subminor) \
34  PETSC_VERSION_EQ(major,minor,subminor)
35 
36 // We used to have workarounds for missing extern "C" in old PETSc
37 // versions. We no longer support PETSc versions so old, but we do
38 // still support libMesh applications old enough to have used these
39 // macros.
40 #define EXTERN_C_FOR_PETSC_BEGIN
41 #define EXTERN_C_FOR_PETSC_END
42 
43 // Petsc include files
44 // Wrapped to avoid triggering our more paranoid warnings
45 #include <libmesh/ignore_warnings.h>
46 #ifdef I
47 # define LIBMESH_SAW_I
48 #endif
49 #include <petsc.h>
50 #ifndef LIBMESH_SAW_I
51 # undef I // Avoid complex.h contamination
52 #endif
53 #include <libmesh/restore_warnings.h>
54 
55 // These macros are still around for backwards compatibility, but we
56 // no longer use them within the library.
57 #define LibMeshVecDestroy(x) VecDestroy(x)
58 #define LibMeshVecScatterDestroy(x) VecScatterDestroy(x)
59 #define LibMeshMatDestroy(x) MatDestroy(x)
60 #define LibMeshISDestroy(x) ISDestroy(x)
61 #define LibMeshKSPDestroy(x) KSPDestroy(x)
62 #define LibMeshSNESDestroy(x) SNESDestroy(x)
63 #define LibMeshPetscViewerDestroy(x) PetscViewerDestroy(x)
64 #define LibMeshPCDestroy(x) PCDestroy(x)
65 
66 // PETSc devs temporarily considered adding VecScatterCreateWithData, but
67 // it was dropped in PETSc-130e142e39 and never made it into any release.
68 // We will keep the ifdef for backwards compatibility in case anyone wrote
69 // code directly using it, but that should be pretty unlikely.
70 #define LibMeshVecScatterCreate(xin,ix,yin,iy,newctx) VecScatterCreate(xin,ix,yin,iy,newctx)
71 #define ISCreateLibMesh(comm,n,idx,mode,is) ISCreateGeneral((comm),(n),(idx),(mode),(is))
72 
73 // As of release 3.8.0, MatGetSubMatrix was renamed to MatCreateSubMatrix.
74 #if PETSC_VERSION_LESS_THAN(3,8,0)
75 # define LibMeshCreateSubMatrix MatGetSubMatrix
76 #else
77 # define LibMeshCreateSubMatrix MatCreateSubMatrix
78 #endif
79 
80 // As of release 3.17.0, SETERRQ was made to use __VA_ARGS__ and
81 // SETERRQ1, SETERRQ2, etc. were deprecated
82 #if PETSC_VERSION_LESS_THAN(3,17,0)
83 # define LIBMESH_SETERRQ1 SETERRQ1
84 # define LIBMESH_SETERRQ2 SETERRQ2
85 # define LIBMESH_SETERRQ3 SETERRQ3
86 #else
87 # define LIBMESH_SETERRQ1 SETERRQ
88 # define LIBMESH_SETERRQ2 SETERRQ
89 # define LIBMESH_SETERRQ3 SETERRQ
90 #endif
91 
92 // As of 3.18, %D is no longer supported in format strings, but the
93 // replacement PetscInt_FMT didn't get added until 3.7.2
94 #if PETSC_VERSION_LESS_THAN(3,8,0)
95 # define LIBMESH_PETSCINT_FMT "D"
96 #else
97 # define LIBMESH_PETSCINT_FMT PetscInt_FMT
98 #endif
99 
100 // As of 3.19, PETSC_NULL is deprecated
101 #if PETSC_VERSION_LESS_THAN(3,19,0)
102 # define LIBMESH_PETSC_NULLPTR PETSC_NULL
103 # define LIBMESH_PETSC_SUCCESS static_cast<PetscErrorCode>(0)
104 #else
105 # define LIBMESH_PETSC_NULLPTR PETSC_NULLPTR
106 # define LIBMESH_PETSC_SUCCESS PETSC_SUCCESS
107 #endif
108 
109 // If we're using quad precision, we need to disambiguate std
110 // operations on PetscScalar
111 
112 #if LIBMESH_DEFAULT_QUADRUPLE_PRECISION
113 # include <boost/multiprecision/float128.hpp>
114 
115 namespace std
116 {
117 
118 // With Boost 1.70 and earlier, we have to define operator<< ourselves
119 // or we can't use it. With Boost 1.71 and later, we *can't* define
120 // operator<< ourselves or it conflicts with theirs.
121 
122 #if BOOST_VERSION < 107000
123 inline
124 std::ostream & operator<< (std::ostream & os, const PetscScalar in)
125 {
126  os << (boost::multiprecision::float128(in));
127  return os;
128 }
129 #endif
130 
131 /*
132 
133 // Whether these are necessary or not depends on gcc version!?
134 
135 #define LIBMESH_PETSCSCALAR_UNARY(funcname) \
136 inline PetscScalar funcname \
137  (const PetscScalar in) \
138 { \
139  return boost::multiprecision::funcname \
140  (boost::multiprecision::float128(in)).backend().value(); \
141 }
142 
143 LIBMESH_PETSCSCALAR_UNARY(sqrt)
144 LIBMESH_PETSCSCALAR_UNARY(exp)
145 LIBMESH_PETSCSCALAR_UNARY(log)
146 LIBMESH_PETSCSCALAR_UNARY(log10)
147 LIBMESH_PETSCSCALAR_UNARY(sin)
148 LIBMESH_PETSCSCALAR_UNARY(cos)
149 LIBMESH_PETSCSCALAR_UNARY(tan)
150 LIBMESH_PETSCSCALAR_UNARY(asin)
151 LIBMESH_PETSCSCALAR_UNARY(acos)
152 LIBMESH_PETSCSCALAR_UNARY(atan)
153 LIBMESH_PETSCSCALAR_UNARY(sinh)
154 LIBMESH_PETSCSCALAR_UNARY(cosh)
155 LIBMESH_PETSCSCALAR_UNARY(tanh)
156 LIBMESH_PETSCSCALAR_UNARY(abs)
157 LIBMESH_PETSCSCALAR_UNARY(fabs)
158 LIBMESH_PETSCSCALAR_UNARY(ceil)
159 LIBMESH_PETSCSCALAR_UNARY(floor)
160 */
161 
162 } // namespace std
163 
164 // Helper functions for boost float128 compatibility
165 namespace libMesh
166 {
167 template <typename T>
168 PetscScalar PS(T val)
169 {
170  return val.backend().value();
171 }
172 
173 template <typename T>
174 PetscScalar * pPS(T * ptr)
175 {
176  return &(ptr->backend().value());
177 }
178 
179 template <typename T>
180 const PetscScalar * pPS(const T * ptr)
181 {
182  return &(ptr->backend().value());
183 }
184 
185 template <typename T>
186 PetscReal * pPR(T * ptr)
187 {
188  return &(ptr->backend().value());
189 }
190 
191 template <typename T>
192 const PetscReal * pPR(const T * ptr)
193 {
194  return &(ptr->backend().value());
195 }
196 } // namespace libMesh
197 
198 #else
199 
200 namespace libMesh
201 {
202 template <typename T>
203 PetscScalar PS(T val)
204 {
205  return val;
206 }
207 
208 template <typename T>
209 PetscScalar * pPS(T * ptr)
210 {
211  return ptr;
212 }
213 
214 template <typename T>
215 const PetscScalar * pPS(const T * ptr)
216 {
217  return ptr;
218 }
219 
220 template <typename T>
221 PetscReal * pPR(T * ptr)
222 {
223  return ptr;
224 }
225 
226 template <typename T>
227 const PetscReal * pPR(const T * ptr)
228 {
229  return ptr;
230 }
231 } // namespace libMesh
232 
233 #endif // LIBMESH_ENABLE_QUADRUPLE_PRECISION
234 
235 #else // LIBMESH_HAVE_PETSC
236 
237 #define PETSC_VERSION_LESS_THAN(major,minor,subminor) 1
238 #define PETSC_VERSION_EQUALS(major,minor,revision) 0
239 
240 #endif // LIBMESH_HAVE_PETSC
241 
242 #define PETSC_RELEASE_LESS_THAN(major, minor, subminor) \
243  PETSC_VERSION_LESS_THAN(major, minor, subminor)
244 #define PETSC_RELEASE_EQUALS(major, minor, subminor) \
245  PETSC_VERSION_EQUALS(major, minor, subminor)
246 
247 #define PETSC_RELEASE_LESS_EQUALS(major, minor, subminor) \
248  (PETSC_RELEASE_LESS_THAN(major, minor, subminor) || PETSC_RELEASE_EQUALS(major, minor, subminor))
249 
250 #define PETSC_RELEASE_GREATER_EQUALS(major, minor, subminor) \
251  (0 == PETSC_RELEASE_LESS_THAN(major, minor, subminor))
252 
253 #define PETSC_RELEASE_GREATER_THAN(major, minor, subminor) \
254  (0 == PETSC_RELEASE_LESS_EQUALS(major, minor, subminor))
255 
256 #endif // LIBMESH_PETSC_MACRO_H
PetscReal * pPR(T *ptr)
Definition: petsc_macro.h:186
std::ostream & operator<<(std::ostream &os, const PetscScalar in)
Definition: petsc_macro.h:124
PetscScalar * pPS(T *ptr)
Definition: petsc_macro.h:174
The libMesh namespace provides an interface to certain functionality in the library.
PetscScalar PS(T val)
Definition: petsc_macro.h:168