GRASS GIS 8 Programmer's Manual  8.4.0dev(2024)-81301443e3
N_pde.h
Go to the documentation of this file.
1 /*****************************************************************************
2  *
3  * MODULE: Grass PDE Numerical Library
4  * AUTHOR(S): Soeren Gebbert, Berlin (GER) Dec 2006
5  * soerengebbert <at> gmx <dot> de
6  *
7  * PURPOSE: This file contains definitions of variables and data types
8  *
9  * COPYRIGHT: (C) 2000 by the GRASS Development Team
10  *
11  * This program is free software under the GNU General Public
12  * License (>=v2). Read the file COPYING that comes with GRASS
13  * for details.
14  *
15  *****************************************************************************/
16 
17 #include <grass/gis.h>
18 #include <grass/raster3d.h>
19 #include <grass/glocale.h>
20 #include <grass/gmath.h>
21 
22 #ifndef _N_PDE_H_
23 #define _N_PDE_H_
24 
25 #define N_NORMAL_LES 0
26 #define N_SPARSE_LES 1
27 /*!
28  * Boundary conditions for cells
29  */
30 #define N_CELL_INACTIVE 0
31 #define N_CELL_ACTIVE 1
32 #define N_CELL_DIRICHLET 2
33 #define N_CELL_TRANSMISSION 3
34 /*!
35  * \brief the maximum number of available cell states (eg: boundary condition,
36  * inactiven active)
37  * */
38 #define N_MAX_CELL_STATE 20
39 
40 #define N_5_POINT_STAR 0
41 #define N_7_POINT_STAR 1
42 #define N_9_POINT_STAR 2
43 #define N_27_POINT_STAR 3
44 
45 #define N_MAXIMUM_NORM 0
46 #define N_EUKLID_NORM 1
47 
48 #define N_ARRAY_SUM 0 /* summ two arrays */
49 #define N_ARRAY_DIF 1 /* calc the difference between two arrays */
50 #define N_ARRAY_MUL 2 /* multiply two arrays */
51 #define N_ARRAY_DIV \
52  3 /* array division, if div with 0 the NULL value is set \
53  */
54 
55 #define N_UPWIND_FULL 0 /*full upwinding stabilization */
56 #define N_UPWIND_EXP 1 /*exponential upwinding stabilization */
57 #define N_UPWIND_WEIGHT 2 /*weighted upwinding stabilization */
58 
59 /* *************************************************************** */
60 /* *************** LINEARE EQUATION SYSTEM PART ****************** */
61 /* *************************************************************** */
62 
63 /*!
64  * \brief The linear equation system (les) structure
65  *
66  * This structure manages the Ax = b system.
67  * It manages regular quadratic matrices or
68  * sparse matrices. The vector b and x are normal one dimensional
69  * memory structures of type double. Also the number of rows
70  * and the matrix type are stored in this structure.
71  * */
72 typedef struct {
73  double *x; /*the value vector */
74  double *b; /*the right side of Ax = b */
75  double **A; /*the normal quadratic matrix */
76  G_math_spvector **Asp; /*the sparse matrix */
77  int rows; /*number of rows */
78  int cols; /*number of cols */
79  int quad; /*is the matrix quadratic (1-quadratic, 0 not) */
80  int type; /*the type of the les, normal == 0, sparse == 1 */
81 } N_les;
82 
83 extern N_les *N_alloc_les_param(int cols, int rows, int type, int param);
84 extern N_les *N_alloc_les(int rows, int type);
85 extern N_les *N_alloc_les_A(int rows, int type);
86 extern N_les *N_alloc_les_Ax(int rows, int type);
87 extern N_les *N_alloc_les_Ax_b(int rows, int type);
88 extern N_les *N_alloc_nquad_les(int cols, int rows, int type);
89 extern N_les *N_alloc_nquad_les_A(int cols, int rows, int type);
90 extern N_les *N_alloc_nquad_les_Ax(int cols, int rows, int type);
91 extern N_les *N_alloc_nquad_les_Ax_b(int cols, int rows, int type);
92 extern void N_print_les(N_les *les);
93 extern void N_free_les(N_les *les);
94 
95 /* *************************************************************** */
96 /* *************** GEOMETRY INFORMATION ************************** */
97 /* *************************************************************** */
98 
99 /*!
100  * \brief Geometric information about the structured grid
101  * */
102 typedef struct {
103  int planimetric; /*If the projection is not planimetric (0), the array
104  calculation is different for each row */
105  double *area; /* the vector of area values for non-planimetric projection
106  for each row */
107  int dim; /* 2 or 3 */
108 
109  double dx;
110  double dy;
111  double dz;
112 
113  double Az;
114 
115  int depths;
116  int rows;
117  int cols;
118 
119 } N_geom_data;
120 
121 extern N_geom_data *N_alloc_geom_data(void);
122 extern void N_free_geom_data(N_geom_data *geodata);
124  N_geom_data *geodata);
125 extern N_geom_data *N_init_geom_data_2d(struct Cell_head *region,
126  N_geom_data *geodata);
127 extern double N_get_geom_data_area_of_cell(N_geom_data *geom, int row);
128 
129 /* *************************************************************** */
130 /* *************** READING RASTER AND VOLUME DATA **************** */
131 /* *************************************************************** */
132 
133 typedef struct {
134  int type; /* which raster type CELL_TYPE, FCELL_TYPE, DCELL_TYPE */
135  int rows, cols;
136  int rows_intern, cols_intern;
137  int offset; /*number of cols/rows offset at each boundary */
138  CELL *cell_array; /*The data is stored in an one dimensional array
139  internally */
140  FCELL *fcell_array; /*The data is stored in an one dimensional array
141  internally */
142  DCELL *dcell_array; /*The data is stored in an one dimensional array
143  internally */
144 } N_array_2d;
145 
146 extern N_array_2d *N_alloc_array_2d(int cols, int rows, int offset, int type);
147 extern void N_free_array_2d(N_array_2d *data_array);
148 extern int N_get_array_2d_type(N_array_2d *array2d);
149 extern void N_get_array_2d_value(N_array_2d *array2d, int col, int row,
150  void *value);
151 extern CELL N_get_array_2d_c_value(N_array_2d *array2d, int col, int row);
152 extern FCELL N_get_array_2d_f_value(N_array_2d *array2d, int col, int row);
153 extern DCELL N_get_array_2d_d_value(N_array_2d *array2d, int col, int row);
154 extern void N_put_array_2d_value(N_array_2d *array2d, int col, int row,
155  char *value);
156 extern void N_put_array_2d_c_value(N_array_2d *array2d, int col, int row,
157  CELL value);
158 extern void N_put_array_2d_f_value(N_array_2d *array2d, int col, int row,
159  FCELL value);
160 extern void N_put_array_2d_d_value(N_array_2d *array2d, int col, int row,
161  DCELL value);
162 extern int N_is_array_2d_value_null(N_array_2d *array2d, int col, int row);
163 extern void N_put_array_2d_value_null(N_array_2d *array2d, int col, int row);
164 extern void N_print_array_2d(N_array_2d *data);
165 extern void N_print_array_2d_info(N_array_2d *data);
166 extern void N_copy_array_2d(N_array_2d *source, N_array_2d *target);
167 extern double N_norm_array_2d(N_array_2d *array1, N_array_2d *array2, int type);
168 extern N_array_2d *N_math_array_2d(N_array_2d *array1, N_array_2d *array2,
169  N_array_2d *result, int type);
171 extern N_array_2d *N_read_rast_to_array_2d(char *name, N_array_2d *array);
172 extern void N_write_array_2d_to_rast(N_array_2d *array, char *name);
173 extern void N_calc_array_2d_stats(N_array_2d *a, double *min, double *max,
174  double *sum, int *nonzero, int withoffset);
175 
176 typedef struct {
177  int type; /* which raster type FCELL_TYPE, DCELL_TYPE */
178  int rows, cols, depths;
179  int rows_intern, cols_intern, depths_intern;
180  int offset; /*number of cols/rows/depths offset at each boundary */
181  float *fcell_array; /*The data is stored in an one dimensional array
182  internally */
183  double *dcell_array; /*The data is stored in an one dimensional array
184  internally */
185 } N_array_3d;
186 
187 extern N_array_3d *N_alloc_array_3d(int cols, int rows, int depths, int offset,
188  int type);
189 extern void N_free_array_3d(N_array_3d *data_array);
190 extern int N_get_array_3d_type(N_array_3d *array3d);
191 extern void N_get_array_3d_value(N_array_3d *array3d, int col, int row,
192  int depth, void *value);
193 extern float N_get_array_3d_f_value(N_array_3d *array3d, int col, int row,
194  int depth);
195 extern double N_get_array_3d_d_value(N_array_3d *array3d, int col, int row,
196  int depth);
197 extern void N_put_array_3d_value(N_array_3d *array3d, int col, int row,
198  int depth, char *value);
199 extern void N_put_array_3d_f_value(N_array_3d *array3d, int col, int row,
200  int depth, float value);
201 extern void N_put_array_3d_d_value(N_array_3d *array3d, int col, int row,
202  int depth, double value);
203 extern int N_is_array_3d_value_null(N_array_3d *array3d, int col, int row,
204  int depth);
205 extern void N_put_array_3d_value_null(N_array_3d *array3d, int col, int row,
206  int depth);
207 extern void N_print_array_3d(N_array_3d *data);
208 extern void N_print_array_3d_info(N_array_3d *data);
209 extern void N_copy_array_3d(N_array_3d *source, N_array_3d *target);
210 extern double N_norm_array_3d(N_array_3d *array1, N_array_3d *array2, int type);
211 extern N_array_3d *N_math_array_3d(N_array_3d *array1, N_array_3d *array2,
212  N_array_3d *result, int type);
215  int mask);
216 extern void N_write_array_3d_to_rast3d(N_array_3d *array, char *name, int mask);
217 extern void N_calc_array_3d_stats(N_array_3d *a, double *min, double *max,
218  double *sum, int *nonzero, int withoffset);
219 
220 /* *************************************************************** */
221 /* *************** MATRIX ASSEMBLING METHODS ********************* */
222 /* *************************************************************** */
223 /*!
224  * \brief Matrix entries for a mass balance 5/7/9 star system
225  *
226  * Matrix entries for the mass balance of a 5 star system
227  *
228  * The entries are center, east, west, north, south and the
229  * right side vector b of Ax = b. This system is typically used in 2d.
230 
231  \verbatim
232  N
233  |
234  W-- C --E
235  |
236  S
237  \endverbatim
238 
239  * Matrix entries for the mass balance of a 7 star system
240  *
241  * The entries are center, east, west, north, south, top, bottom and the
242  * right side vector b of Ax = b. This system is typically used in 3d.
243 
244  \verbatim
245  T N
246  |/
247  W-- C --E
248  /|
249  S B
250  \endverbatim
251 
252  * Matrix entries for the mass balance of a 9 star system
253  *
254  * The entries are center, east, west, north, south, north-east, south-east,
255  * north-wast, south-west and the
256  * right side vector b of Ax = b. This system is typically used in 2d.
257 
258  \verbatim
259  NW N NE
260  \ | /
261  W-- C --E
262  / | \
263  SW S SE
264  \endverbatim
265 
266  * Matrix entries for the mass balance of a 27 star system
267  *
268  * The entries are center, east, west, north, south, north-east, south-east,
269  * north-wast, south-west, same for top and bottom and the
270  * right side vector b of Ax = b. This system is typically used in 2d.
271 
272  \verbatim
273  top:
274  NW_T N_Z NE_T
275  \ | /
276  W_T-- T --E_T
277  / | \
278  SW_T S_T SE_T
279 
280  center:
281  NW N NE
282  \ | /
283  W-- C --E
284  / | \
285  SW S SE
286 
287  bottom:
288  NW_B N_B NE_B
289  \ | /
290  W_B-- B --E_B
291  / | \
292  SW_B S_B SE_B
293  \endverbatim
294 
295  */
296 typedef struct {
297  int type;
298  int count;
299  double C, W, E, N, S, NE, NW, SE, SW, V;
300  /*top part */
301  double T, W_T, E_T, N_T, S_T, NE_T, NW_T, SE_T, SW_T;
302  /*bottom part */
303  double B, W_B, E_B, N_B, S_B, NE_B, NW_B, SE_B, SW_B;
304 } N_data_star;
305 
306 /*!
307  * \brief callback structure for 3d matrix assembling
308  * */
309 typedef struct {
310  N_data_star *(*callback)(void *, N_geom_data *, int, int, int);
312 
313 /*!
314  * \brief callback structure for 2d matrix assembling
315  * */
316 typedef struct {
317  N_data_star *(*callback)(void *, N_geom_data *, int, int);
319 
320 extern void N_set_les_callback_3d_func(
321  N_les_callback_3d *data,
322  N_data_star *(*callback_func_3d)(void *, N_geom_data *, int, int, int));
323 extern void N_set_les_callback_2d_func(
324  N_les_callback_2d *data,
325  N_data_star *(*callback_func_2d)(void *, N_geom_data *, int, int));
328 extern N_data_star *N_alloc_5star(void);
329 extern N_data_star *N_alloc_7star(void);
330 extern N_data_star *N_alloc_9star(void);
331 extern N_data_star *N_alloc_27star(void);
332 extern N_data_star *N_create_5star(double C, double W, double E, double N,
333  double S, double V);
334 extern N_data_star *N_create_7star(double C, double W, double E, double N,
335  double S, double T, double B, double V);
336 extern N_data_star *N_create_9star(double C, double W, double E, double N,
337  double S, double NW, double SW, double NE,
338  double SE, double V);
339 extern N_data_star *
340 N_create_27star(double C, double W, double E, double N, double S, double NW,
341  double SW, double NE, double SE, double T, double W_T,
342  double E_T, double N_T, double S_T, double NW_T, double SW_T,
343  double NE_T, double SE_T, double B, double W_B, double E_B,
344  double N_B, double S_B, double NW_B, double SW_B, double NE_B,
345  double SE_B, double V);
346 extern N_data_star *N_callback_template_3d(void *data, N_geom_data *geom,
347  int col, int row, int depth);
348 extern N_data_star *N_callback_template_2d(void *data, N_geom_data *geom,
349  int col, int row);
350 extern N_les *N_assemble_les_3d(int les_type, N_geom_data *geom,
351  N_array_3d *status, N_array_3d *start_val,
352  void *data, N_les_callback_3d *callback);
353 extern N_les *N_assemble_les_3d_active(int les_type, N_geom_data *geom,
354  N_array_3d *status,
355  N_array_3d *start_val, void *data,
356  N_les_callback_3d *callback);
357 extern N_les *N_assemble_les_3d_dirichlet(int les_type, N_geom_data *geom,
358  N_array_3d *status,
359  N_array_3d *start_val, void *data,
360  N_les_callback_3d *callback);
361 extern N_les *N_assemble_les_3d_param(int les_type, N_geom_data *geom,
362  N_array_3d *status, N_array_3d *start_val,
363  void *data, N_les_callback_3d *callback,
364  int cell_type);
365 extern N_les *N_assemble_les_2d(int les_type, N_geom_data *geom,
366  N_array_2d *status, N_array_2d *start_val,
367  void *data, N_les_callback_2d *callback);
368 extern N_les *N_assemble_les_2d_active(int les_type, N_geom_data *geom,
369  N_array_2d *status,
370  N_array_2d *start_val, void *data,
371  N_les_callback_2d *callback);
372 extern N_les *N_assemble_les_2d_dirichlet(int les_type, N_geom_data *geom,
373  N_array_2d *status,
374  N_array_2d *start_val, void *data,
375  N_les_callback_2d *callback);
376 extern N_les *N_assemble_les_2d_param(int les_type, N_geom_data *geom,
377  N_array_2d *status, N_array_2d *start_val,
378  void *data, N_les_callback_2d *callback,
379  int cell_Type);
380 extern int N_les_pivot_create(N_les *les);
382  N_array_2d *status, N_array_2d *start_val);
384  N_array_3d *status, N_array_3d *start_val);
385 
386 /* *************************************************************** */
387 /* *************** GPDE STANDARD OPTIONS ************************* */
388 /* *************************************************************** */
389 
390 /*! \brief Standard options of the gpde library
391  * */
392 typedef enum {
393  N_OPT_SOLVER_SYMM, /*! solver for symmetric, positive definite linear
394  equation systems */
395  N_OPT_SOLVER_UNSYMM, /*! solver for unsymmetric linear equation systems */
396  N_OPT_MAX_ITERATIONS, /*! Maximum number of iteration used to solver the
397  linear equation system */
398  N_OPT_ITERATION_ERROR, /*! Error break criteria for the iterative solver
399  (jacobi, sor, cg or bicgstab) */
400  N_OPT_SOR_VALUE, /*! The relaxation parameter used by the jacobi and sor
401  solver for speedup or stabilizing */
402  N_OPT_CALC_TIME /*! The calculation time in seconds */
403 } N_STD_OPT;
404 
405 extern struct Option *N_define_standard_option(int opt);
406 
407 /* *************************************************************** */
408 /* *************** GPDE MATHEMATICAL TOOLS *********************** */
409 /* *************************************************************** */
410 
411 extern double N_calc_arith_mean(double a, double b);
412 extern double N_calc_arith_mean_n(double *a, int size);
413 extern double N_calc_geom_mean(double a, double b);
414 extern double N_calc_geom_mean_n(double *a, int size);
415 extern double N_calc_harmonic_mean(double a, double b);
416 extern double N_calc_harmonic_mean_n(double *a, int size);
417 extern double N_calc_quad_mean(double a, double b);
418 extern double N_calc_quad_mean_n(double *a, int size);
419 
420 /* *************************************************************** */
421 /* *************** UPWIND STABILIZATION ALGORITHMS *************** */
422 /* *************************************************************** */
423 
424 extern double N_full_upwinding(double sprod, double distance, double D);
425 extern double N_exp_upwinding(double sprod, double distance, double D);
426 
427 /* *************************************************************** */
428 /* *************** METHODS FOR GRADIENT CALCULATION ************** */
429 /* *************************************************************** */
430 /*!
431  \verbatim
432 
433  ______________
434  | | | |
435  | | | |
436  |----|-NC-|----|
437  | | | |
438  | WC EC |
439  | | | |
440  |----|-SC-|----|
441  | | | |
442  |____|____|____|
443 
444 
445  | /
446  TC NC
447  |/
448  --WC-----EC--
449  /|
450  SC BC
451  / |
452 
453  \endverbatim
454 
455  */
456 
457 /*! \brief Gradient between the cells in X and Y direction */
458 typedef struct {
459 
460  double NC, SC, WC, EC;
461 
462 } N_gradient_2d;
463 
464 /*! \brief Gradient between the cells in X, Y and Z direction */
465 typedef struct {
466 
467  double NC, SC, WC, EC, TC, BC;
468 
469 } N_gradient_3d;
470 
471 /*!
472  \verbatim
473 
474  Gradient in X direction between the cell neighbours
475  ____ ____ ____
476  | | | |
477  | NWN NEN |
478  |____|____|____|
479  | | | |
480  | WN EN |
481  |____|____|____|
482  | | | |
483  | SWS SES |
484  |____|____|____|
485 
486  Gradient in Y direction between the cell neighbours
487  ______________
488  | | | |
489  | | | |
490  |NWW-|-NC-|-NEE|
491  | | | |
492  | | | |
493  |SWW-|-SC-|-SEE|
494  | | | |
495  |____|____|____|
496 
497  Gradient in Z direction between the cell neighbours
498  /______________/
499  /| | | |
500  | NWZ| NZ | NEZ|
501  |____|____|____|
502  /| | | |
503  | WZ | CZ | EZ |
504  |____|____|____|
505  /| | | |
506  | SWZ| SZ | SEZ|
507  |____|____|____|
508  /____/____/____/
509 
510 
511  \endverbatim
512  */
513 
514 /*! \brief Gradient between the cell neighbours in X direction */
515 typedef struct {
516 
517  double NWN, NEN, WC, EC, SWS, SES;
518 
520 
521 /*! \brief Gradient between the cell neighbours in Y direction */
522 typedef struct {
523 
524  double NWW, NEE, NC, SC, SWW, SEE;
525 
527 
528 /*! \brief Gradient between the cell neighbours in Z direction */
529 typedef struct {
530 
531  double NWZ, NZ, NEZ, WZ, CZ, EZ, SWZ, SZ, SEZ;
532 
534 
535 /*! \brief Gradient between the cell neighbours in X and Y direction */
536 typedef struct {
537 
540 
542 
543 /*! \brief Gradient between the cell neighbours in X, Y and Z direction */
544 typedef struct {
545 
546  N_gradient_neighbours_x *xt; /*top values */
547  N_gradient_neighbours_x *xc; /*center values */
548  N_gradient_neighbours_x *xb; /*bottom values */
549 
550  N_gradient_neighbours_y *yt; /*top values */
551  N_gradient_neighbours_y *yc; /*center values */
552  N_gradient_neighbours_y *yb; /*bottom values */
553 
554  N_gradient_neighbours_z *zt; /*top-center values */
555  N_gradient_neighbours_z *zb; /*bottom-center values */
556 
558 
559 /*! Two dimensional gradient field */
560 typedef struct {
561 
562  N_array_2d *x_array;
563  N_array_2d *y_array;
564  int cols, rows;
565  double min, max, mean, sum;
566  int nonull;
567 
569 
570 /*! Three dimensional gradient field */
571 typedef struct {
572 
573  N_array_3d *x_array;
574  N_array_3d *y_array;
575  N_array_3d *z_array;
576  int cols, rows, depths;
577  double min, max, mean, sum;
578  int nonull;
579 
581 
582 extern N_gradient_2d *N_alloc_gradient_2d(void);
583 extern void N_free_gradient_2d(N_gradient_2d *grad);
584 extern N_gradient_2d *N_create_gradient_2d(double NC, double SC, double WC,
585  double EC);
586 extern int N_copy_gradient_2d(N_gradient_2d *source, N_gradient_2d *target);
588  N_gradient_2d *gradient, int col,
589  int row);
590 extern N_gradient_3d *N_alloc_gradient_3d(void);
591 extern void N_free_gradient_3d(N_gradient_3d *grad);
592 extern N_gradient_3d *N_create_gradient_3d(double NC, double SC, double WC,
593  double EC, double TC, double BC);
594 extern int N_copy_gradient_3d(N_gradient_3d *source, N_gradient_3d *target);
596  N_gradient_3d *gradient, int col,
597  int row, int depth);
601 N_create_gradient_neighbours_x(double NWN, double NEN, double WC, double EC,
602  double SWS, double SES);
604  N_gradient_neighbours_x *target);
608 N_create_gradient_neighbours_y(double NWW, double NEE, double NC, double SC,
609  double SWW, double SEE);
611  N_gradient_neighbours_y *target);
615 N_create_gradient_neighbours_z(double NWZ, double NZ, double NEZ, double WZ,
616  double CZ, double EZ, double SWZ, double SZ,
617  double SEZ);
619  N_gradient_neighbours_z *target);
626  N_gradient_neighbours_2d *target);
629  N_gradient_neighbours_2d *gradient, int col,
630  int row);
639  N_gradient_neighbours_3d *target);
642 extern N_gradient_field_2d *N_alloc_gradient_field_2d(int cols, int rows);
645  N_gradient_field_2d *target);
646 extern N_gradient_field_2d *
648  N_array_2d *weight_y, N_geom_data *geom,
649  N_gradient_field_2d *gradfield);
651  N_array_2d *x_comp,
652  N_array_2d *y_comp);
655 extern N_gradient_field_3d *N_alloc_gradient_field_3d(int cols, int rows,
656  int depths);
659  N_gradient_field_3d *target);
660 extern N_gradient_field_3d *
662  N_array_3d *weight_y, N_array_3d *weight_z,
663  N_geom_data *geom, N_gradient_field_3d *gradfield);
665  N_array_3d *x_comp,
666  N_array_3d *y_comp,
667  N_array_3d *z_comp);
668 
669 #endif
N_data_star * N_alloc_27star(void)
allocate a 27 point star data structure
double N_norm_array_2d(N_array_2d *array1, N_array_2d *array2, int type)
Calculate the norm of the two input arrays.
N_les * N_assemble_les_2d_active(int les_type, N_geom_data *geom, N_array_2d *status, N_array_2d *start_val, void *data, N_les_callback_2d *callback)
Assemble a linear equation system (les) based on 2d location data (raster) and active cells.
void N_print_les(N_les *les)
prints the linear equation system to stdout
Definition: n_les.c:257
int N_copy_gradient_neighbours_2d(N_gradient_neighbours_2d *source, N_gradient_neighbours_2d *target)
copy a N_gradient_neighbours_2d structure
Definition: n_gradient.c:650
N_gradient_neighbours_y * N_alloc_gradient_neighbours_y(void)
Allocate a N_gradient_neighbours_y structure.
Definition: n_gradient.c:381
void N_compute_gradient_field_components_3d(N_gradient_field_3d *field, N_array_3d *x_comp, N_array_3d *y_comp, N_array_3d *z_comp)
Calculate the x, y and z vector components from a gradient field for each cell and store them in the ...
N_gradient_neighbours_3d * N_alloc_gradient_neighbours_3d(void)
Allocate a N_gradient_neighbours_3d structure.
Definition: n_gradient.c:752
void N_put_array_3d_f_value(N_array_3d *array3d, int col, int row, int depth, float value)
This function writes a float value to the N_array_3d data at position col, row, depth.
Definition: n_arrays.c:1121
CELL N_get_array_2d_c_value(N_array_2d *array2d, int col, int row)
Returns the value of type CELL at position col, row.
Definition: n_arrays.c:314
int N_copy_gradient_3d(N_gradient_3d *source, N_gradient_3d *target)
copy a N_gradient_3d structure
Definition: n_gradient.c:211
N_gradient_2d * N_alloc_gradient_2d(void)
Allocate a N_gradient_2d structure.
Definition: n_gradient.c:26
void N_print_array_3d_info(N_array_3d *data)
Write the info of the array to stdout.
Definition: n_arrays.c:1170
double N_get_geom_data_area_of_cell(N_geom_data *geom, int row)
Get the areay size in square meter of one cell (x*y) at row.
Definition: n_geom.c:196
N_data_star * N_alloc_9star(void)
allocate a 9 point star data structure
N_les * N_assemble_les_2d_param(int les_type, N_geom_data *geom, N_array_2d *status, N_array_2d *start_val, void *data, N_les_callback_2d *callback, int cell_Type)
Assemble a linear equation system (les) based on 2d location data (raster)
DCELL N_get_array_2d_d_value(N_array_2d *array2d, int col, int row)
Returns the value of type DCELL at position col, row.
Definition: n_arrays.c:380
N_gradient_3d * N_get_gradient_3d(N_gradient_field_3d *field, N_gradient_3d *gradient, int col, int row, int depth)
Return a N_gradient_3d structure calculated from the input gradient field at position [depth][row][co...
Definition: n_gradient.c:246
void N_put_array_3d_value_null(N_array_3d *array3d, int col, int row, int depth)
This function writes a null value to the N_array_3d data at position col, row, depth.
Definition: n_arrays.c:1060
N_STD_OPT
Standard options of the gpde library.
Definition: N_pde.h:391
@ N_OPT_MAX_ITERATIONS
Definition: N_pde.h:395
@ N_OPT_SOLVER_UNSYMM
Definition: N_pde.h:394
@ N_OPT_SOLVER_SYMM
Definition: N_pde.h:392
@ N_OPT_ITERATION_ERROR
Definition: N_pde.h:397
@ N_OPT_CALC_TIME
Definition: N_pde.h:401
@ N_OPT_SOR_VALUE
Definition: N_pde.h:399
N_gradient_2d * N_create_gradient_2d(double NC, double SC, double WC, double EC)
allocate and initialize a N_gradient_2d structure
Definition: n_gradient.c:59
N_gradient_neighbours_2d * N_alloc_gradient_neighbours_2d(void)
Allocate a N_gradient_neighbours_2d structure.
Definition: n_gradient.c:577
void N_put_array_3d_d_value(N_array_3d *array3d, int col, int row, int depth, double value)
Writes a double value to the N_array_3d struct at position col, row, depth.
Definition: n_arrays.c:1148
void N_compute_gradient_field_components_2d(N_gradient_field_2d *field, N_array_2d *x_comp, N_array_2d *y_comp)
Calculate the x and y vector components from a gradient field for each cell and stores them in the pr...
void N_calc_gradient_field_2d_stats(N_gradient_field_2d *field)
Calculate basic statistics of a gradient field.
void N_put_array_3d_value(N_array_3d *array3d, int col, int row, int depth, char *value)
This function writes a value to the N_array_3d data at position col, row, depth.
Definition: n_arrays.c:1010
void N_free_gradient_field_3d(N_gradient_field_3d *field)
Free's a N_gradient_neighbours_3d structure.
Definition: n_gradient.c:1019
void N_write_array_3d_to_rast3d(N_array_3d *array, char *name, int mask)
Write a N_array_3d struct to a volume map.
Definition: n_arrays_io.c:377
N_les * N_assemble_les_3d(int les_type, N_geom_data *geom, N_array_3d *status, N_array_3d *start_val, void *data, N_les_callback_3d *callback)
Assemble a linear equation system (les) based on 3d location data (g3d) active cells.
int N_convert_array_3d_null_to_zero(N_array_3d *a)
Convert all null values to zero values.
double N_norm_array_3d(N_array_3d *array1, N_array_3d *array2, int type)
Calculate the norm of the two input arrays.
N_array_3d * N_read_rast3d_to_array_3d(char *name, N_array_3d *array, int mask)
Read a volume map into a N_array_3d structure.
Definition: n_arrays_io.c:249
double N_exp_upwinding(double sprod, double distance, double D)
exponential upwinding stabilization algorithm
Definition: n_upwind.c:63
int N_convert_array_2d_null_to_zero(N_array_2d *a)
Convert all null values to zero values.
N_gradient_neighbours_2d * N_create_gradient_neighbours_2d(N_gradient_neighbours_x *x, N_gradient_neighbours_y *y)
Allocate and initialize a N_gradient_neighbours_2d structure.
Definition: n_gradient.c:618
int N_copy_gradient_field_3d(N_gradient_field_3d *source, N_gradient_field_3d *target)
Copy N_gradient_field_3d structure from source to target.
Definition: n_gradient.c:1041
int N_copy_gradient_neighbours_x(N_gradient_neighbours_x *source, N_gradient_neighbours_x *target)
copy a N_gradient_neighbours_x structure
Definition: n_gradient.c:355
int N_les_integrate_dirichlet_3d(N_les *les, N_geom_data *geom, N_array_3d *status, N_array_3d *start_val)
Integrate Dirichlet or Transmission boundary conditions into the les (3d)
int N_get_array_3d_type(N_array_3d *array3d)
Return the data type of the N_array_3d.
Definition: n_arrays.c:803
FCELL N_get_array_2d_f_value(N_array_2d *array2d, int col, int row)
Returns the value of type FCELL at position col, row.
Definition: n_arrays.c:347
N_data_star * N_callback_template_2d(void *data, N_geom_data *geom, int col, int row)
N_array_2d * N_read_rast_to_array_2d(char *name, N_array_2d *array)
Read a raster map into a N_array_2d structure.
Definition: n_arrays_io.c:44
N_geom_data * N_init_geom_data_2d(struct Cell_head *region, N_geom_data *geodata)
Initiate a pde geometry data structure with a 2d region.
Definition: n_geom.c:117
N_gradient_field_3d * N_compute_gradient_field_3d(N_array_3d *pot, N_array_3d *weight_x, N_array_3d *weight_y, N_array_3d *weight_z, N_geom_data *geom, N_gradient_field_3d *gradfield)
This function computes the gradient based on the input N_array_3d pot (that means potential),...
void N_free_gradient_neighbours_2d(N_gradient_neighbours_2d *grad)
Free's a N_gradient_neighbours_2d structure.
Definition: n_gradient.c:596
N_les * N_alloc_nquad_les(int cols, int rows, int type)
Allocate memory for a (not) quadratic linear equation system which includes the Matrix A,...
Definition: n_les.c:34
N_les * N_alloc_nquad_les_A(int cols, int rows, int type)
Allocate memory for a (not) quadratic linear equation system which includes the Matrix A.
Definition: n_les.c:68
void N_print_gradient_field_2d_info(N_gradient_field_2d *field)
Print gradient field information to stdout.
Definition: n_gradient.c:961
double N_calc_geom_mean_n(double *a, int size)
Calculate the geometrical mean of the values in vector a of size n.
Definition: n_tools.c:93
N_gradient_neighbours_x * N_alloc_gradient_neighbours_x(void)
Allocate a N_gradient_neighbours_x structure.
Definition: n_gradient.c:289
N_les_callback_3d * N_alloc_les_callback_3d(void)
Allocate the structure holding the callback function.
N_gradient_neighbours_y * N_create_gradient_neighbours_y(double NWW, double NEE, double NC, double SC, double SWW, double SEE)
Allocate and initialize a N_gradient_neighbours_y structure.
Definition: n_gradient.c:418
int N_is_array_2d_value_null(N_array_2d *array2d, int col, int row)
Returns 1 if the value of N_array_2d struct at position col, row is of type null, otherwise 0.
Definition: n_arrays.c:231
double N_calc_harmonic_mean_n(double *a, int size)
Calculate the harmonical mean of the values in vector a of size n.
Definition: n_tools.c:136
void N_print_gradient_field_3d_info(N_gradient_field_3d *field)
Print gradient field information to stdout.
Definition: n_gradient.c:1062
void N_free_gradient_neighbours_z(N_gradient_neighbours_z *grad)
Free's a N_gradient_neighbours_z structure.
Definition: n_gradient.c:489
int N_copy_gradient_neighbours_z(N_gradient_neighbours_z *source, N_gradient_neighbours_z *target)
copy a N_gradient_neighbours_z structure
Definition: n_gradient.c:547
void N_free_array_3d(N_array_3d *data_array)
Release the memory of a N_array_3d.
Definition: n_arrays.c:774
N_data_star * N_create_9star(double C, double W, double E, double N, double S, double NW, double SW, double NE, double SE, double V)
allocate and initialize a 9 point star data structure
N_array_3d * N_alloc_array_3d(int cols, int rows, int depths, int offset, int type)
Allocate memory for a N_array_3d data structure.
Definition: n_arrays.c:719
N_data_star * N_callback_template_3d(void *data, N_geom_data *geom, int col, int row, int depth)
N_les * N_alloc_les_Ax(int rows, int type)
Allocate memory for a quadratic linear equation system which includes the Matrix A and vector x.
Definition: n_les.c:117
N_les * N_alloc_les(int rows, int type)
Allocate memory for a quadratic linear equation system which includes the Matrix A,...
Definition: n_les.c:101
int N_les_integrate_dirichlet_2d(N_les *les, N_geom_data *geom, N_array_2d *status, N_array_2d *start_val)
Integrate Dirichlet or Transmission boundary conditions into the les (2s)
void N_calc_array_3d_stats(N_array_3d *a, double *min, double *max, double *sum, int *nonzero, int withoffset)
Calculate basic statistics of the N_array_3d struct.
N_geom_data * N_init_geom_data_3d(RASTER3D_Region *region3d, N_geom_data *geodata)
Initiate a pde geometry data structure with a 3d region.
Definition: n_geom.c:73
void N_print_array_2d(N_array_2d *data)
Write info and content of the N_array_2d struct to stdout.
Definition: n_arrays.c:628
int N_copy_gradient_neighbours_3d(N_gradient_neighbours_3d *source, N_gradient_neighbours_3d *target)
copy a N_gradient_neighbours_3d structure
Definition: n_gradient.c:852
N_les * N_assemble_les_2d(int les_type, N_geom_data *geom, N_array_2d *status, N_array_2d *start_val, void *data, N_les_callback_2d *callback)
Assemble a linear equation system (les) based on 2d location data (raster) and active cells.
N_les_callback_2d * N_alloc_les_callback_2d(void)
Allocate the structure holding the callback function.
void N_free_geom_data(N_geom_data *geodata)
Release memory of a pde geometry data structure.
Definition: n_geom.c:49
double N_calc_quad_mean_n(double *a, int size)
Calculate the quadratic mean of the values in vector a of size n.
Definition: n_tools.c:184
void N_get_array_2d_value(N_array_2d *array2d, int col, int row, void *value)
Write the value of the N_array_2d struct at position col, row to value.
Definition: n_arrays.c:182
N_gradient_field_2d * N_alloc_gradient_field_2d(int cols, int rows)
Allocate a N_gradient_field_2d.
Definition: n_gradient.c:896
double N_calc_geom_mean(double a, double b)
Calculate the geometrical mean of values a and b.
Definition: n_tools.c:73
double N_calc_arith_mean_n(double *a, int size)
Calculate the arithmetic mean of the values in vector a of size n.
Definition: n_tools.c:51
void N_copy_array_3d(N_array_3d *source, N_array_3d *target)
Copy the source N_array_3d struct to the target N_array_3d struct.
void N_print_array_2d_info(N_array_2d *data)
This function writes the data info of the array data to stdout.
Definition: n_arrays.c:603
double N_calc_arith_mean(double a, double b)
Calculate the arithmetic mean of values a and b.
Definition: n_tools.c:31
N_gradient_neighbours_2d * N_get_gradient_neighbours_2d(N_gradient_field_2d *field, N_gradient_neighbours_2d *gradient, int col, int row)
Return a N_gradient_neighbours_2d structure calculated from the input gradient field at position [row...
Definition: n_gradient.c:689
N_gradient_3d * N_create_gradient_3d(double NC, double SC, double WC, double EC, double TC, double BC)
allocate and initialize a N_gradient_3d structure
Definition: n_gradient.c:184
N_array_2d * N_math_array_2d(N_array_2d *array1, N_array_2d *array2, N_array_2d *result, int type)
Perform calculations with two input arrays, the result is written to a third array.
void N_copy_array_2d(N_array_2d *source, N_array_2d *target)
Copy the source N_array_2d struct to the target N_array_2d struct.
Definition: n_arrays_calc.c:43
N_gradient_3d * N_alloc_gradient_3d(void)
Allocate a N_gradient_3d structure.
Definition: n_gradient.c:149
N_gradient_2d * N_get_gradient_2d(N_gradient_field_2d *field, N_gradient_2d *gradient, int col, int row)
Return a N_gradient_2d structure calculated from the input gradient field at position [row][col].
Definition: n_gradient.c:114
N_array_2d * N_alloc_array_2d(int cols, int rows, int offset, int type)
Allocate memory for a N_array_2d data structure.
Definition: n_arrays.c:75
N_les * N_assemble_les_3d_active(int les_type, N_geom_data *geom, N_array_3d *status, N_array_3d *start_val, void *data, N_les_callback_3d *callback)
Assemble a linear equation system (les) based on 3d location data (g3d) active cells.
N_les * N_alloc_nquad_les_Ax(int cols, int rows, int type)
Allocate memory for a (not) quadratic linear equation system which includes the Matrix A and vector x...
Definition: n_les.c:51
int N_les_pivot_create(N_les *les)
double N_full_upwinding(double sprod, double distance, double D)
full upwinding stabilization algorithm
Definition: n_upwind.c:32
void N_free_gradient_2d(N_gradient_2d *grad)
Free's a N_gradient_2d structure.
Definition: n_gradient.c:41
int N_copy_gradient_neighbours_y(N_gradient_neighbours_y *source, N_gradient_neighbours_y *target)
copy a N_gradient_neighbours_y structure
Definition: n_gradient.c:447
void N_free_les(N_les *les)
Release the memory of the linear equation system.
Definition: n_les.c:307
int N_is_array_3d_value_null(N_array_3d *array3d, int col, int row, int depth)
This function returns 1 if value of N_array_3d data at position col, row, depth is of type null,...
Definition: n_arrays.c:873
N_array_3d * N_math_array_3d(N_array_3d *array1, N_array_3d *array2, N_array_3d *result, int type)
Perform calculations with two input arrays, the result is written to a third array.
N_les * N_assemble_les_3d_param(int les_type, N_geom_data *geom, N_array_3d *status, N_array_3d *start_val, void *data, N_les_callback_3d *callback, int cell_type)
Assemble a linear equation system (les) based on 3d location data (g3d)
N_gradient_neighbours_z * N_alloc_gradient_neighbours_z(void)
Allocate a N_gradient_neighbours_z structure.
Definition: n_gradient.c:473
void N_put_array_2d_d_value(N_array_2d *array2d, int col, int row, DCELL value)
Writes a DCELL value to the N_array_2d struct at position col, row.
Definition: n_arrays.c:576
N_data_star * N_alloc_7star(void)
allocate a 7 point star data structure
int N_copy_gradient_2d(N_gradient_2d *source, N_gradient_2d *target)
copy a N_gradient_2d structure
Definition: n_gradient.c:83
int N_get_array_2d_type(N_array_2d *array2d)
Return the data type of the N_array_2d struct.
Definition: n_arrays.c:164
N_geom_data * N_alloc_geom_data(void)
Allocate the pde geometry data structure and return a pointer to the new allocated structure.
Definition: n_geom.c:29
void N_free_gradient_3d(N_gradient_3d *grad)
Free's a N_gradient_3d structure.
Definition: n_gradient.c:164
void N_write_array_2d_to_rast(N_array_2d *array, char *name)
Write a N_array_2d struct to a raster map.
Definition: n_arrays_io.c:171
void N_put_array_2d_value_null(N_array_2d *array2d, int col, int row)
Writes the null value to the N_array_2d struct at position col, row.
Definition: n_arrays.c:459
N_data_star * N_create_27star(double C, double W, double E, double N, double S, double NW, double SW, double NE, double SE, double T, double W_T, double E_T, double N_T, double S_T, double NW_T, double SW_T, double NE_T, double SE_T, double B, double W_B, double E_B, double N_B, double S_B, double NW_B, double SW_B, double NE_B, double SE_B, double V)
allocate and initialize a 27 point star data structure
void N_calc_array_2d_stats(N_array_2d *a, double *min, double *max, double *sum, int *nonzero, int withoffset)
Calculate basic statistics of the N_array_2d struct.
void N_free_gradient_neighbours_y(N_gradient_neighbours_y *grad)
Free's a N_gradient_neighbours_y structure.
Definition: n_gradient.c:397
N_gradient_field_2d * N_compute_gradient_field_2d(N_array_2d *pot, N_array_2d *weight_x, N_array_2d *weight_y, N_geom_data *geom, N_gradient_field_2d *gradfield)
This function computes the gradient based on the input N_array_2d pot (potential),...
N_les * N_assemble_les_2d_dirichlet(int les_type, N_geom_data *geom, N_array_2d *status, N_array_2d *start_val, void *data, N_les_callback_2d *callback)
Assemble a linear equation system (les) based on 2d location data (raster) and active and dirichlet c...
struct Option * N_define_standard_option(int opt)
Create standardised Option structure related to the gpde library.
void N_get_array_3d_value(N_array_3d *array3d, int col, int row, int depth, void *value)
This function writes the value of N_array_3d data at position col, row, depth to the variable value.
Definition: n_arrays.c:822
void N_calc_gradient_field_3d_stats(N_gradient_field_3d *field)
Calculate basic statistics of a gradient field.
N_gradient_neighbours_x * N_create_gradient_neighbours_x(double NWN, double NEN, double WC, double EC, double SWS, double SES)
Allocate and initialize a N_gradient_neighbours_x structure.
Definition: n_gradient.c:326
N_gradient_field_3d * N_alloc_gradient_field_3d(int cols, int rows, int depths)
Allocate a N_gradient_field_3d.
Definition: n_gradient.c:993
N_les * N_alloc_les_Ax_b(int rows, int type)
Allocate memory for a quadratic linear equation system which includes the Matrix A,...
Definition: n_les.c:149
double N_calc_quad_mean(double a, double b)
Calculate the quadratic mean of values a and b.
Definition: n_tools.c:164
void N_put_array_2d_f_value(N_array_2d *array2d, int col, int row, FCELL value)
Writes a FCELL value to the N_array_2d struct at position col, row.
Definition: n_arrays.c:546
float N_get_array_3d_f_value(N_array_3d *array3d, int col, int row, int depth)
This function returns the value of type float at position col, row, depth.
Definition: n_arrays.c:948
N_les * N_alloc_les_A(int rows, int type)
Allocate memory for a quadratic linear equation system which includes the Matrix A.
Definition: n_les.c:133
N_data_star * N_alloc_5star(void)
allocate a 5 point star data structure
double N_calc_harmonic_mean(double a, double b)
Calculate the harmonical mean of values a and b.
Definition: n_tools.c:115
double N_get_array_3d_d_value(N_array_3d *array3d, int col, int row, int depth)
This function returns the value of type float at position col, row, depth.
Definition: n_arrays.c:979
N_les * N_assemble_les_3d_dirichlet(int les_type, N_geom_data *geom, N_array_3d *status, N_array_3d *start_val, void *data, N_les_callback_3d *callback)
Assemble a linear equation system (les) based on 3d location data (g3d) active and dirichlet cells.
void N_free_array_2d(N_array_2d *data_array)
Release the memory of a N_array_2d structure.
Definition: n_arrays.c:132
void N_set_les_callback_3d_func(N_les_callback_3d *data, N_data_star *(*callback_func_3d)(void *, N_geom_data *, int, int, int))
Set the callback function which is called while assembling the les in 3d.
void N_put_array_2d_c_value(N_array_2d *array2d, int col, int row, CELL value)
Writes a CELL value to the N_array_2d struct at position col, row.
Definition: n_arrays.c:516
void N_print_array_3d(N_array_3d *data)
Write info and content of the array data to stdout.
Definition: n_arrays.c:1196
void N_set_les_callback_2d_func(N_les_callback_2d *data, N_data_star *(*callback_func_2d)(void *, N_geom_data *, int, int))
Set the callback function which is called while assembling the les in 2d.
void N_free_gradient_field_2d(N_gradient_field_2d *field)
Free's a N_gradient_neighbours_2d structure.
Definition: n_gradient.c:920
N_data_star * N_create_7star(double C, double W, double E, double N, double S, double T, double B, double V)
allocate and initialize a 7 point star data structure
N_gradient_neighbours_3d * N_create_gradient_neighbours_3d(N_gradient_neighbours_x *xt, N_gradient_neighbours_x *xc, N_gradient_neighbours_x *xb, N_gradient_neighbours_y *yt, N_gradient_neighbours_y *yc, N_gradient_neighbours_y *yb, N_gradient_neighbours_z *zt, N_gradient_neighbours_z *zb)
Allocate and initialize a N_gradient_neighbours_3d structure.
Definition: n_gradient.c:806
N_data_star * N_create_5star(double C, double W, double E, double N, double S, double V)
allocate and initialize a 5 point star data structure
void N_free_gradient_neighbours_x(N_gradient_neighbours_x *grad)
Free's a N_gradient_neighbours_x structure.
Definition: n_gradient.c:305
N_gradient_neighbours_z * N_create_gradient_neighbours_z(double NWZ, double NZ, double NEZ, double WZ, double CZ, double EZ, double SWZ, double SZ, double SEZ)
Allocate and initialize a N_gradient_neighbours_z structure.
Definition: n_gradient.c:513
N_les * N_alloc_les_param(int cols, int rows, int type, int param)
Allocate memory for a quadratic or not quadratic linear equation system.
Definition: n_les.c:182
int N_copy_gradient_field_2d(N_gradient_field_2d *source, N_gradient_field_2d *target)
Copy N_gradient_field_2d structure from source to target.
Definition: n_gradient.c:941
N_les * N_alloc_nquad_les_Ax_b(int cols, int rows, int type)
Allocate memory for a (not) quadratic linear equation system which includes the Matrix A,...
Definition: n_les.c:85
void N_free_gradient_neighbours_3d(N_gradient_neighbours_3d *grad)
Free's a N_gradient_neighbours_3d structure.
Definition: n_gradient.c:777
void N_put_array_2d_value(N_array_2d *array2d, int col, int row, char *value)
Writes a value to the N_array_2d struct at position col, row.
Definition: n_arrays.c:412
#define SE
Definition: dataquad.h:31
#define SW
Definition: dataquad.h:30
#define NE
Definition: dataquad.h:29
#define NW
Definition: dataquad.h:28
#define min(x, y)
Definition: draw2.c:29
#define max(x, y)
Definition: draw2.c:30
#define N
Definition: e_intersect.c:926
#define D
Definition: gis/intersect.c:72
float FCELL
Definition: gis.h:627
double DCELL
Definition: gis.h:626
int CELL
Definition: gis.h:625
float mean(IClass_statistics *statistics, int band)
Helper function for computing mean.
int count
const char * name
Definition: named_colr.c:6
#define W
Definition: ogsf.h:143
double b
Definition: r_raster.c:39
2D/3D raster map header (used also for region)
Definition: gis.h:437
The row vector of the sparse matrix.
Definition: gmath.h:54
Matrix entries for a mass balance 5/7/9 star system.
Definition: N_pde.h:295
Geometric information about the structured grid.
Definition: N_pde.h:101
Gradient between the cells in X and Y direction.
Definition: N_pde.h:457
Gradient between the cells in X, Y and Z direction.
Definition: N_pde.h:464
Gradient between the cell neighbours in X and Y direction.
Definition: N_pde.h:535
Gradient between the cell neighbours in X, Y and Z direction.
Definition: N_pde.h:543
Gradient between the cell neighbours in X direction.
Definition: N_pde.h:514
Gradient between the cell neighbours in Y direction.
Definition: N_pde.h:521
Gradient between the cell neighbours in Z direction.
Definition: N_pde.h:528
callback structure for 2d matrix assembling
Definition: N_pde.h:315
callback structure for 3d matrix assembling
Definition: N_pde.h:308
The linear equation system (les) structure.
Definition: N_pde.h:71
Structure that stores option information.
Definition: gis.h:554
#define x