1 #ifndef _PYGSL_SOLVER_H_
2 #define _PYGSL_SOLVER_H_ 1
3 #include <pygsl/intern.h>
4 #include <pygsl/block_helpers.h>
5 /* Not directly needed here, but provides a lot of convienience functions */
6 #include <pygsl/error_helpers.h>
7 #include <pygsl/string_helpers.h>
8 #include <gsl/gsl_math.h>
9 #include <setjmp.h>
10 
11 #undef __BEGIN_DECLS
12 #undef __END_DECLS
13 #ifdef __cplusplus
14 # define __BEGIN_DECLS extern "C" {
15 # define __END_DECLS }
16 #else
17 # define __BEGIN_DECLS /* empty */
18 # define __END_DECLS /* empty */
19 #endif
20 
21 __BEGIN_DECLS
22 
23 /*
24  * Many functions are "just" accessor methods. These different methods are
25  * listed here.
26  *
27  * Convention: they all end with "m_t" short for method type. This emphasises
28  * that this type is pointer to a method.
29  */
30 typedef int (*int_m_t)(void *);
31 typedef size_t (*size_t_m_t)(void *);
32 typedef void (*void_m_t)(void *);
33 typedef void *(*void_a_t)(const void *);
34 typedef void *(*void_an_t)(const void *, size_t n);
35 typedef void *(*void_anp_t)(const void *, size_t n, size_t p);
36 typedef const char * (*name_m_t)(void *);
37 typedef double  (*double_m_t)(void *);
38 typedef gsl_vector *(*ret_vec)(void *);
39 typedef int (*set_m_t)(void *, void *, const gsl_vector *);
40 typedef int (*set_m_d_t) (void *, gsl_function *, double);
41 typedef int (*set_m_ddd_t) (void *, gsl_function *, double, double, double);
42 typedef int (*int_f_vd_t)(const gsl_vector *, double);
43 typedef int (*int_f_vvdd_t)(const gsl_vector *, const gsl_vector *, double, double);
44 
45 struct _PyGSLSolverObject;
46 
47 /*
48  * GSL Methods which are implemented as C functions. The specifiy pointer to
49  * the solver struct are passed as (void *) pointers.
50  */
51 struct _GSLMethods{
52      /* Method to be called to free the GSL solver */
53      void_m_t free;
54      /*
55       * Some solvers provide a restart method. If not available set it to NULL.
56       */
57      void_m_t restart;
58      /* returns a string with the name of the solver */
59      name_m_t name;
60      /* takes one more step towards the solution */
61      int_m_t iterate;
62 };
63 
64 struct _SolverStatic{
65      struct _GSLMethods cmethods;
66      /* How many callbacks will be used? */
67      int n_cbs;
68      /* Additional methods not provided by the basis solver. */
69      PyMethodDef *pymethods;
70      /* Describes the type of the solver. e. g. F-Minimizer */
71      const char * type_name;
72 };
73 
74 struct pygsl_array_cache{
75      double * data;
76      PyArrayObject * ref;
77 };
78 #define PyGSL_SOLVER_NCBS_MAX 4
79 #define PyGSL_SOLVER_PB_ND_MAX 2
80 #define PyGSL_SOLVER_N_ARRAYS 10
81 
82 struct _PyGSLSolverObject{
83      PyObject_HEAD
84      /*
85       *	Some solvers do not propagate errors. Here I use longjmp if an error
86       * is raised by the evaluator.
87       */
88      jmp_buf buffer;
89      /*
90       * Array generation is a vital part of the callback process. To increase
91       * the calculation speed, I want to store them here, thus I can reuse them
92       * if not stored by the user for later evaluation.
93       */
94      struct pygsl_array_cache *cache;
95      /*
96       * The Python callback methods.
97       */
98      PyObject* cbs[PyGSL_SOLVER_NCBS_MAX];
99      /*
100       * Additional arguments passed to the callbacks
101       */
102      PyObject* args;
103      /*
104       * The solver itself.
105       */
106      void * solver;
107      /*
108       * The space needed to store the variables for the C functions ....
109       */
110      void * c_sys;
111      /*
112       * The dimensionality of the problem. Typically one or two numbers ...
113       */
114      int problem_dimensions[PyGSL_SOLVER_PB_ND_MAX];
115      /*
116       * The methods of the solver.
117       */
118      const struct _SolverStatic* mstatic;
119      /* before one can iterate the solver the method set() must be called!*/
120      int set_called;
121      /* Used as a flag if the jmp_buf is set */
122      int isset;
123 };
124 typedef struct _PyGSLSolverObject PyGSL_solver;
125 
126 
127 typedef struct{
128      const void * type;
129      void* alloc;
130      const struct _SolverStatic* mstatic;
131 } solver_alloc_struct;
132 
133 
134 
135 #ifndef _PyGSL_SOLVER_API_MODULE
136 #define PyGSL_SOLVER_API_EXTERN extern
137 #else
138 #define PyGSL_SOLVER_API_EXTERN static
139 #endif
140 
141 /*
142  *  Initalises a solver
143  *
144  *  The internal structure is set up and the
145  *  nd: number of dimensions
146  *     0 ... zero dimensional e.g. minimize, root
147  *     1 ... one  dimensional e.g. multimin
148  *     2 ... two  dimensional e.g. multifit_nlin
149  *
150  *     3 ... only initalise the structure
151  *
152  */
153 PyGSL_SOLVER_API_EXTERN  PyObject *
154 PyGSL_solver_dn_init(PyObject *self, PyObject *args, const solver_alloc_struct * alloc, int nd);
155 
156 /*
157  * Accessor methods.
158  *
159  * These methods allow to access parameters of the C structure using the access
160  * methods _m_t func
161  */
162 PyGSL_SOLVER_API_EXTERN PyObject*
163 PyGSL_solver_ret_double(PyGSL_solver *self, PyObject *args, double_m_t func);
164 
165 PyGSL_SOLVER_API_EXTERN PyObject*
166 PyGSL_solver_ret_int(PyGSL_solver *self, PyObject *args, int_m_t func);
167 
168 PyGSL_SOLVER_API_EXTERN PyObject*
169 PyGSL_solver_ret_size_t(PyGSL_solver *self, PyObject *args, size_t_m_t func);
170 
171 PyGSL_SOLVER_API_EXTERN PyObject*
172 PyGSL_solver_ret_vec(PyGSL_solver *self, PyObject *args,  ret_vec func);
173 
174 /*
175  * evaluates a C function taking an vector and a double as input and returning a status.
176  */
177 
178 #define PyGSL_CALLABLE_CHECK(ob, name) \
179 (PyCallable_Check(ob) ? GSL_SUCCESS : PyGSL_Callable_Check(ob, name))
180 
181 PyGSL_SOLVER_API_EXTERN int
182 PyGSL_Callable_Check(PyObject *f, const char * myname);
183 
184 PyGSL_SOLVER_API_EXTERN int
185 PyGSL_function_wrap_On_O(const gsl_vector * x, PyObject *callback,
186 			 PyObject *arguments, double *result1,
187 			 gsl_vector *result2, int n, const char * c_func_name);
188 
189 PyGSL_SOLVER_API_EXTERN int
190 PyGSL_function_wrap_OnOn_On(const gsl_vector *x, const gsl_vector *v, gsl_vector *hv, PyObject *callback,
191 			    PyObject *arguments, int n, const char *c_func_name);
192 
193 PyGSL_SOLVER_API_EXTERN int
194 PyGSL_function_wrap_Op_On_Opn(const gsl_vector *x, gsl_vector *f1, gsl_matrix *f2, PyObject *callback,
195 			    PyObject *arguments, int n, int p, const char *c_func_name);
196 
197 PyGSL_SOLVER_API_EXTERN int
198 PyGSL_function_wrap_Op_On(const gsl_vector * x, gsl_vector *f, PyObject *callback,
199 			  PyObject * arguments, int n, int p, const char *c_func_name);
200 
201 PyGSL_SOLVER_API_EXTERN int
202 PyGSL_function_wrap_Op_Opn(const gsl_vector * x, gsl_matrix *f, PyObject *callback,
203 			   PyObject *arguments, int n, int p, const char * c_func_name);
204 
205 /*
206  * evaluates a C function taking an vector and a double as input and returning a status.
207  */
208 PyGSL_SOLVER_API_EXTERN PyObject*
209 PyGSL_solver_vd_i(PyObject * self, PyObject *args, int_f_vd_t func);
210 
211 PyGSL_SOLVER_API_EXTERN PyObject *
212 PyGSL_solver_vvdd_i(PyObject * self, PyObject * args, int_f_vvdd_t func);
213 
214 struct pygsl_solver_n_set{
215      int is_fdf;
216      void *c_sys;
217      set_m_t set;
218 };
219 
220 PyGSL_SOLVER_API_EXTERN PyObject *
221 PyGSL_solver_n_set(PyGSL_solver *self, PyObject *pyargs, PyObject *kw,
222 		   const struct pygsl_solver_n_set * info);
223 
224 
225 /*
226  *
227  */
228 PyGSL_SOLVER_API_EXTERN int
229 PyGSL_solver_func_set(PyGSL_solver *self, PyObject *args, PyObject *f,
230 		       PyObject *df, PyObject *fdf);
231 
232 PyGSL_SOLVER_API_EXTERN PyObject*
233 PyGSL_solver_set_f(PyGSL_solver *self, PyObject *pyargs, PyObject *kw,
234 		   void *fptr, int isfdf);
235 
236 
237 #define _GET(name, cast, func) \
238 PyObject * PyGSL_ ## name(PyGSL_solver *self, PyObject *args) \
239 { \
240     return PyGSL_solver_ ## func(self, args, (cast) gsl_ ## name); \
241 }
242 
243 #define GETDOUBLE(name) _GET(name, double_m_t, ret_double)
244 #define GETSIZET(name)  _GET(name, size_t_m_t, ret_size_t)
245 #define GETINT(name)    _GET(name, int_m_t,    ret_int)
246 #define GETVEC(name)    _GET(name, ret_vec,    ret_vec)
247 
248 
249 /*
250  * Get or set a double ....
251  */
252 enum PyGSL_GETSET_typemode {
253      PyGSL_MODE_DOUBLE = 0,
254      PyGSL_MODE_INT,
255      PyGSL_MODE_SIZE_T
256 };
257 
258 PyGSL_SOLVER_API_EXTERN PyObject *
259 PyGSL_solver_GetSet(PyObject *self, PyObject *args, void * address, enum PyGSL_GETSET_typemode mode);
260 
261 #ifndef _PyGSL_SOLVER_API_MODULE
262 
263 #define PyGSL_function_wrap_Op_On \
264 (* (int (*)(const gsl_vector *, gsl_vector *, PyObject *, PyObject *, int, int, const char *))\
265  PyGSL_API[PyGSL_function_wrap_Op_On_NUM])
266 
267 #define PyGSL_function_wrap_On_O \
268 (* (int (*)(const gsl_vector *, PyObject *, PyObject *, double *, gsl_vector *, int, const char *))\
269  PyGSL_API[PyGSL_function_wrap_On_O_NUM])
270 
271 #define PyGSL_function_wrap_OnOn_On \
272 (* (int (*)(const gsl_vector *, const gsl_vector *, gsl_vector *, PyObject *,  PyObject *,int, const char *))\
273  PyGSL_API[PyGSL_function_wrap_OnOn_On_NUM])
274 
275 #define PyGSL_function_wrap_Op_On_Opn \
276 (* (int (*)(const gsl_vector *, gsl_vector *, gsl_matrix *, PyObject *, PyObject *, int, int, const char *))\
277  PyGSL_API[PyGSL_function_wrap_Op_On_Opn_NUM])
278 
279 #define PyGSL_function_wrap_Op_Opn \
280 (* (int (*)(const gsl_vector *, gsl_matrix *, PyObject *, PyObject *, int, int, const char *))\
281  PyGSL_API[PyGSL_function_wrap_Op_Opn_NUM])
282 
283 
284 #define PyGSL_solver_ret_int \
285 (*(PyObject *  (*) (PyGSL_solver *, PyObject *, int_m_t func))    PyGSL_API[PyGSL_solver_ret_int_NUM])
286 
287 #define PyGSL_solver_ret_double \
288 (*(PyObject *  (*) (PyGSL_solver *, PyObject *, double_m_t func)) PyGSL_API[PyGSL_solver_ret_double_NUM])
289 
290 #define PyGSL_solver_ret_size_t \
291 (*(PyObject *  (*) (PyGSL_solver *, PyObject *, size_t_m_t func)) PyGSL_API[PyGSL_solver_ret_size_t_NUM])
292 
293 #define PyGSL_solver_ret_vec \
294 (*(PyObject *  (*) (PyGSL_solver *, PyObject *, ret_vec func))    PyGSL_API[PyGSL_solver_ret_vec_NUM])
295 
296 #define  PyGSL_solver_dn_init \
297 (*(PyObject * (*)(PyObject *, PyObject *, const solver_alloc_struct *, int))PyGSL_API[PyGSL_solver_dn_init_NUM])
298 
299 #define PyGSL_Callable_Check \
300 (*(int (*)(PyObject *, const char *)) PyGSL_API[PyGSL_Callable_Check_NUM])
301 
302 #define PyGSL_solver_vd_i \
303 (*(PyObject * (*) (PyObject *, PyObject *, int_f_vd_t)) PyGSL_API[PyGSL_solver_vd_i_NUM])
304 
305 #define PyGSL_solver_vvdd_i \
306 (*(PyObject * (*) (PyObject *, PyObject *, int_f_vvdd_t)) PyGSL_API[PyGSL_solver_vvdd_i_NUM])
307 
308 #define PyGSL_solver_func_set \
309 (*(int (*)(PyGSL_solver *, PyObject *, PyObject *, PyObject *, PyObject *)) PyGSL_API[PyGSL_solver_func_set_NUM])
310 
311 #define PyGSL_solver_n_set \
312 (*(PyObject * (*)(PyGSL_solver *, PyObject *, PyObject *,  const struct pygsl_solver_n_set *)) PyGSL_API[PyGSL_solver_n_set_NUM])
313 
314 #define PyGSL_solver_set_f \
315 (* (PyObject * (*)(PyGSL_solver *, PyObject *, PyObject *, void *, int )) PyGSL_API[PyGSL_solver_set_f_NUM])
316 
317 #define PyGSL_solver_GetSet \
318 (* (PyObject * (*) (PyObject *, PyObject *, void *, enum PyGSL_GETSET_typemode mode)) PyGSL_API[PyGSL_solver_getset_NUM])
319 
320 #define PyGSL_solver_check(ob) ((Py_TYPE(ob)) == (PyGSL_API[PyGSL_solver_type_NUM]))
321 
322 #define import_pygsl_solver() \
323 { \
324    init_pygsl(); \
325    if (PyImport_ImportModule("pygsl.testing.solver") != NULL) { \
326           ;\
327    } else { \
328         fprintf(stderr, "failed to import pygsl solver!!\n"); \
329    } \
330 }
331 
332 #else  /* _PyGSL_API_MODULE */
333 #define PyGSL_solver_check(ob) (Py_TYPE(ob) == &PyGSL_solver_pytype)
334 #endif /* _PyGSL_API_MODULE */
335 
336 __END_DECLS
337 
338 #endif /* _PYGSL_SOLVER_H_ */
339