GRASS Programmer's Manual  6.4.4(2014)-r
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Macros Pages
N_pde.h
Go to the documentation of this file.
1 
2 /*****************************************************************************
3 *
4 * MODULE: Grass PDE Numerical Library
5 * AUTHOR(S): Soeren Gebbert, Berlin (GER) Dec 2006
6 * soerengebbert <at> gmx <dot> de
7 *
8 * PURPOSE: This file contains definitions of variables and data types
9 *
10 * COPYRIGHT: (C) 2000 by the GRASS Development Team
11 *
12 * This program is free software under the GNU General Public
13 * License (>=v2). Read the file COPYING that comes with GRASS
14 * for details.
15 *
16 *****************************************************************************/
17 
18 #include <grass/gis.h>
19 #include <grass/G3d.h>
20 #include <grass/glocale.h>
21 
22 #ifndef _N_PDE_H_
23 #define _N_PDE_H_
24 
25 /*solver names */
26 #define N_SOLVER_DIRECT_GAUSS "gauss"
27 #define N_SOLVER_DIRECT_LU "lu"
28 #define N_SOLVER_DIRECT_CHOLESKY "cholesky"
29 #define N_SOLVER_ITERATIVE_JACOBI "jacobi"
30 #define N_SOLVER_ITERATIVE_SOR "sor"
31 #define N_SOLVER_ITERATIVE_CG "cg"
32 #define N_SOLVER_ITERATIVE_PCG "pcg"
33 #define N_SOLVER_ITERATIVE_BICGSTAB "bicgstab"
34 
35 /*preconditioner */
36 #define N_DIAGONAL_PRECONDITION 1
37 #define N_ROWSCALE_ABSSUMNORM_PRECONDITION 2
38 #define N_ROWSCALE_EUKLIDNORM_PRECONDITION 3
39 #define N_ROWSCALE_MAXNORM_PRECONDITION 4
40 
41 #define N_NORMAL_LES 0
42 #define N_SPARSE_LES 1
43 
44 #define N_CELL_INACTIVE 0
45 #define N_CELL_ACTIVE 1
46 #define N_CELL_DIRICHLET 2
47 #define N_CELL_TRANSMISSION 3
48 
51 #define N_MAX_CELL_STATE 20
52 
53 #define N_5_POINT_STAR 0
54 #define N_7_POINT_STAR 1
55 #define N_9_POINT_STAR 2
56 #define N_27_POINT_STAR 3
57 
58 #define N_MAXIMUM_NORM 0
59 #define N_EUKLID_NORM 1
60 
61 #define N_ARRAY_SUM 0 /* summ two arrays */
62 #define N_ARRAY_DIF 1 /* calc the difference between two arrays */
63 #define N_ARRAY_MUL 2 /* multiply two arrays */
64 #define N_ARRAY_DIV 3 /* array division, if div with 0 the NULL value is set */
65 
66 #define N_UPWIND_FULL 0 /*full upwinding stabilization */
67 #define N_UPWIND_EXP 1 /*exponential upwinding stabilization */
68 #define N_UPWIND_WEIGHT 2 /*weighted upwinding stabilization */
69 
70 
71 
72 /* *************************************************************** */
73 /* *************** LINEARE EQUATION SYSTEM PART ****************** */
74 /* *************************************************************** */
75 
79 typedef struct
80 {
81  int cols; /*Number of entries */
82  double *values; /*The non null values of the row */
83  int *index; /*the index number */
84 } N_spvector;
85 
86 
96 typedef struct
97 {
98  double *x; /*the value vector */
99  double *b; /*the right side of Ax = b */
100  double **A; /*the normal quadratic matrix */
101  N_spvector **Asp; /*the sparse matrix */
102  int rows; /*number of rows */
103  int cols; /*number of cols */
104  int quad; /*is the matrix quadratic (1-quadratic, 0 not) */
105  int type; /*the type of the les, normal == 0, sparse == 1 */
106 } N_les;
107 
108 extern N_spvector *N_alloc_spvector(int cols);
109 extern N_les *N_alloc_les_param(int cols, int rows, int type, int param);
110 extern N_les *N_alloc_les(int rows, int type);
111 extern N_les *N_alloc_les_A(int rows, int type);
112 extern N_les *N_alloc_les_Ax(int rows, int type);
113 extern N_les *N_alloc_les_Ax_b(int rows, int type);
114 extern N_les *N_alloc_nquad_les(int cols, int rows, int type);
115 extern N_les *N_alloc_nquad_les_A(int cols, int rows, int type);
116 extern N_les *N_alloc_nquad_les_Ax(int cols, int rows, int type);
117 extern N_les *N_alloc_nquad_les_Ax_b(int cols, int rows, int type);
118 extern void N_print_les(N_les * les);
119 extern int N_add_spvector_to_les(N_les * les, N_spvector * vector, int row);
120 extern void N_free_spvector(N_spvector * vector);
121 extern void N_free_les(N_les * les);
122 
123 /* *************************************************************** */
124 /* *************** GEOMETRY INFORMATION ************************** */
125 /* *************************************************************** */
126 
130 typedef struct
131 {
132  int planimetric; /*If the projection is not planimetric (0), the array calculation is different for each row */
133  double *area; /* the vector of area values for non-planimetric projection for each row */
134  int dim; /* 2 or 3 */
135 
136  double dx;
137  double dy;
138  double dz;
139 
140  double Az;
141 
142  int depths;
143  int rows;
144  int cols;
145 
146 } N_geom_data;
147 
148 extern N_geom_data *N_alloc_geom_data(void);
149 extern void N_free_geom_data(N_geom_data * geodata);
150 extern N_geom_data *N_init_geom_data_3d(G3D_Region * region3d,
151  N_geom_data * geodata);
152 extern N_geom_data *N_init_geom_data_2d(struct Cell_head *region,
153  N_geom_data * geodata);
154 extern double N_get_geom_data_area_of_cell(N_geom_data * geom, int row);
155 
156 
157 /* *************************************************************** */
158 /* *************** LINEARE EQUATION SOLVER PART ****************** */
159 /* *************************************************************** */
160 extern int N_solver_gauss(N_les * les);
161 extern int N_solver_lu(N_les * les);
162 extern int N_solver_cholesky(N_les * les);
163 extern int N_solver_jacobi(N_les * L, int maxit, double sor, double error);
164 extern int N_solver_SOR(N_les * L, int maxit, double sor, double error);
165 extern int N_solver_cg(N_les * les, int maxit, double error);
166 extern int N_solver_pcg(N_les * les, int maxit, double error, int prec);
167 extern int N_solver_bicgstab(N_les * les, int maxit, double error);
168 extern void N_matrix_vector_product(N_les * les, double *source,
169  double *result);
170 extern void N_sparse_matrix_vector_product(N_les * les, double *source,
171  double *result);
172 extern N_les *N_create_diag_precond_matrix(N_les * les, int prec);
173 
174 /* *************************************************************** */
175 /* *************** READING RASTER AND VOLUME DATA **************** */
176 /* *************************************************************** */
177 
178 typedef struct
179 {
180  int type; /* which raster type CELL_TYPE, FCELL_TYPE, DCELL_TYPE */
181  int rows, cols;
182  int rows_intern, cols_intern;
183  int offset; /*number of cols/rows offset at each boundary */
184  CELL *cell_array; /*The data is stored in an one dimensional array internally */
185  FCELL *fcell_array; /*The data is stored in an one dimensional array internally */
186  DCELL *dcell_array; /*The data is stored in an one dimensional array internally */
187 } N_array_2d;
188 
189 extern N_array_2d *N_alloc_array_2d(int cols, int rows, int offset, int type);
190 extern void N_free_array_2d(N_array_2d * data_array);
191 extern int N_get_array_2d_type(N_array_2d * array2d);
192 extern void N_get_array_2d_value(N_array_2d * array2d, int col, int row,
193  void *value);
194 extern CELL N_get_array_2d_c_value(N_array_2d * array2d, int col, int row);
195 extern FCELL N_get_array_2d_f_value(N_array_2d * array2d, int col, int row);
196 extern DCELL N_get_array_2d_d_value(N_array_2d * array2d, int col, int row);
197 extern void N_put_array_2d_value(N_array_2d * array2d, int col, int row,
198  char *value);
199 extern void N_put_array_2d_c_value(N_array_2d * array2d, int col, int row,
200  CELL value);
201 extern void N_put_array_2d_f_value(N_array_2d * array2d, int col, int row,
202  FCELL value);
203 extern void N_put_array_2d_d_value(N_array_2d * array2d, int col, int row,
204  DCELL value);
205 extern int N_is_array_2d_value_null(N_array_2d * array2d, int col, int row);
206 extern void N_put_array_2d_value_null(N_array_2d * array2d, int col, int row);
207 extern void N_print_array_2d(N_array_2d * data);
208 extern void N_print_array_2d_info(N_array_2d * data);
209 extern void N_copy_array_2d(N_array_2d * source, N_array_2d * target);
210 extern double N_norm_array_2d(N_array_2d * array1, N_array_2d * array2,
211  int type);
212 extern N_array_2d *N_math_array_2d(N_array_2d * array1, N_array_2d * array2,
213  N_array_2d * result, int type);
215 extern N_array_2d *N_read_rast_to_array_2d(char *name, N_array_2d * array);
216 extern void N_write_array_2d_to_rast(N_array_2d * array, char *name);
217 extern void N_calc_array_2d_stats(N_array_2d * a, double *min, double *max,
218  double *sum, int *nonzero, int withoffset);
219 
220 typedef struct
221 {
222  int type; /* which raster type FCELL_TYPE, DCELL_TYPE */
223  int rows, cols, depths;
224  int rows_intern, cols_intern, depths_intern;
225  int offset; /*number of cols/rows/depths offset at each boundary */
226  float *fcell_array; /*The data is stored in an one dimensional array internally */
227  double *dcell_array; /*The data is stored in an one dimensional array internally */
228 } N_array_3d;
229 
230 extern N_array_3d *N_alloc_array_3d(int cols, int rows, int depths,
231  int offset, int type);
232 extern void N_free_array_3d(N_array_3d * data_array);
233 extern int N_get_array_3d_type(N_array_3d * array3d);
234 extern void N_get_array_3d_value(N_array_3d * array3d, int col, int row,
235  int depth, void *value);
236 extern float N_get_array_3d_f_value(N_array_3d * array3d, int col, int row,
237  int depth);
238 extern double N_get_array_3d_d_value(N_array_3d * array3d, int col, int row,
239  int depth);
240 extern void N_put_array_3d_value(N_array_3d * array3d, int col, int row,
241  int depth, char *value);
242 extern void N_put_array_3d_f_value(N_array_3d * array3d, int col, int row,
243  int depth, float value);
244 extern void N_put_array_3d_d_value(N_array_3d * array3d, int col, int row,
245  int depth, double value);
246 extern int N_is_array_3d_value_null(N_array_3d * array3d, int col, int row,
247  int depth);
248 extern void N_put_array_3d_value_null(N_array_3d * array3d, int col, int row,
249  int depth);
250 extern void N_print_array_3d(N_array_3d * data);
251 extern void N_print_array_3d_info(N_array_3d * data);
252 extern void N_copy_array_3d(N_array_3d * source, N_array_3d * target);
253 extern double N_norm_array_3d(N_array_3d * array1, N_array_3d * array2,
254  int type);
255 extern N_array_3d *N_math_array_3d(N_array_3d * array1, N_array_3d * array2,
256  N_array_3d * result, int type);
258 extern N_array_3d *N_read_rast3d_to_array_3d(char *name, N_array_3d * array,
259  int mask);
260 extern void N_write_array_3d_to_rast3d(N_array_3d * array, char *name,
261  int mask);
262 extern void N_calc_array_3d_stats(N_array_3d * a, double *min, double *max,
263  double *sum, int *nonzero, int withoffset);
264 
265 /* *************************************************************** */
266 /* *************** MATRIX ASSEMBLING METHODS ********************* */
267 /* *************************************************************** */
341 typedef struct
342 {
343  int type;
344  int count;
345  double C, W, E, N, S, NE, NW, SE, SW, V;
346  /*top part */
347  double T, W_T, E_T, N_T, S_T, NE_T, NW_T, SE_T, SW_T;
348  /*bottom part */
349  double B, W_B, E_B, N_B, S_B, NE_B, NW_B, SE_B, SW_B;
350 } N_data_star;
351 
355 typedef struct
356 {
357  N_data_star *(*callback) ();
359 
363 typedef struct
364 {
365  N_data_star *(*callback) ();
367 
368 
370  N_data_star * (*callback_func_3d) ());
372  N_data_star * (*callback_func_2d) ());
375 extern N_data_star *N_alloc_5star(void);
376 extern N_data_star *N_alloc_7star(void);
377 extern N_data_star *N_alloc_9star(void);
378 extern N_data_star *N_alloc_27star(void);
379 extern N_data_star *N_create_5star(double C, double W, double E, double N,
380  double S, double V);
381 extern N_data_star *N_create_7star(double C, double W, double E, double N,
382  double S, double T, double B, double V);
383 extern N_data_star *N_create_9star(double C, double W, double E, double N,
384  double S, double NW, double SW, double NE,
385  double SE, double V);
386 extern N_data_star *N_create_27star(double C, double W, double E, double N,
387  double S, double NW, double SW, double NE,
388  double SE, double T, double W_T,
389  double E_T, double N_T, double S_T,
390  double NW_T, double SW_T, double NE_T,
391  double SE_T, double B, double W_B,
392  double E_B, double N_B, double S_B,
393  double NW_B, double SW_B, double NE_B,
394  double SE_B, double V);
395 
396 extern N_data_star *N_callback_template_3d(void *data, N_geom_data * geom,
397  int col, int row, int depth);
398 extern N_data_star *N_callback_template_2d(void *data, N_geom_data * geom,
399  int col, int row);
400 extern N_les *N_assemble_les_3d(int les_type, N_geom_data * geom,
401  N_array_3d * status, N_array_3d * start_val,
402  void *data, N_les_callback_3d * callback);
403 extern N_les *N_assemble_les_3d_active(int les_type, N_geom_data * geom,
404  N_array_3d * status,
405  N_array_3d * start_val, void *data,
406  N_les_callback_3d * callback);
407 extern N_les *N_assemble_les_3d_dirichlet(int les_type, N_geom_data * geom,
408  N_array_3d * status,
409  N_array_3d * start_val, void *data,
410  N_les_callback_3d * callback);
411 extern N_les *N_assemble_les_3d_param(int les_type, N_geom_data * geom,
412  N_array_3d * status,
413  N_array_3d * start_val, void *data,
414  N_les_callback_3d * callback,
415  int cell_type);
416 extern N_les *N_assemble_les_2d(int les_type, N_geom_data * geom,
417  N_array_2d * status, N_array_2d * start_val,
418  void *data, N_les_callback_2d * callback);
419 extern N_les *N_assemble_les_2d_active(int les_type, N_geom_data * geom,
420  N_array_2d * status,
421  N_array_2d * start_val, void *data,
422  N_les_callback_2d * callback);
423 extern N_les *N_assemble_les_2d_dirichlet(int les_type, N_geom_data * geom,
424  N_array_2d * status,
425  N_array_2d * start_val, void *data,
426  N_les_callback_2d * callback);
427 extern N_les *N_assemble_les_2d_param(int les_type, N_geom_data * geom,
428  N_array_2d * status,
429  N_array_2d * start_val, void *data,
430  N_les_callback_2d * callback,
431  int cell_Type);
432 
433 extern int N_les_pivot_create(N_les * les);
435  N_array_2d * status, N_array_2d * start_val);
437  N_array_3d * status, N_array_3d * start_val);
438 
439 /* *************************************************************** */
440 /* *************** GPDE STANDARD OPTIONS ************************* */
441 /* *************************************************************** */
442 
445 typedef enum
446 {
453 } N_STD_OPT;
454 
455 extern struct Option *N_define_standard_option(int opt);
456 
457 /* *************************************************************** */
458 /* *************** GPDE MATHEMATICAL TOOLS *********************** */
459 /* *************************************************************** */
460 
461 extern double N_calc_arith_mean(double a, double b);
462 extern double N_calc_arith_mean_n(double *a, int size);
463 extern double N_calc_geom_mean(double a, double b);
464 extern double N_calc_geom_mean_n(double *a, int size);
465 extern double N_calc_harmonic_mean(double a, double b);
466 extern double N_calc_harmonic_mean_n(double *a, int size);
467 extern double N_calc_quad_mean(double a, double b);
468 extern double N_calc_quad_mean_n(double *a, int size);
469 
470 /* *************************************************************** */
471 /* *************** UPWIND STABILIZATION ALGORITHMS *************** */
472 /* *************************************************************** */
473 
474 extern double N_full_upwinding(double sprod, double distance, double D);
475 extern double N_exp_upwinding(double sprod, double distance, double D);
476 
477 
478 /* *************************************************************** */
479 /* *************** METHODS FOR GRADIENT CALCULATION ************** */
480 /* *************************************************************** */
509 typedef struct
510 {
511 
512  double NC, SC, WC, EC;
513 
514 } N_gradient_2d;
515 
517 typedef struct
518 {
519 
520  double NC, SC, WC, EC, TC, BC;
521 
522 } N_gradient_3d;
523 
524 
569 typedef struct
570 {
571 
572  double NWN, NEN, WC, EC, SWS, SES;
573 
575 
577 typedef struct
578 {
579 
580  double NWW, NEE, NC, SC, SWW, SEE;
581 
583 
585 typedef struct
586 {
587 
588  double NWZ, NZ, NEZ, WZ, CZ, EZ, SWZ, SZ, SEZ;
589 
591 
593 typedef struct
594 {
595 
598 
600 
601 
603 typedef struct
604 {
605 
606  N_gradient_neighbours_x *xt; /*top values */
607  N_gradient_neighbours_x *xc; /*center values */
608  N_gradient_neighbours_x *xb; /*bottom values */
609 
610  N_gradient_neighbours_y *yt; /*top values */
611  N_gradient_neighbours_y *yc; /*center values */
612  N_gradient_neighbours_y *yb; /*bottom values */
613 
614  N_gradient_neighbours_z *zt; /*top-center values */
615  N_gradient_neighbours_z *zb; /*bottom-center values */
616 
618 
619 
621 typedef struct
622 {
623 
626  int cols, rows;
627  double min, max, mean, sum;
628  int nonull;
629 
631 
633 typedef struct
634 {
635 
639  int cols, rows, depths;
640  double min, max, mean, sum;
641  int nonull;
642 
644 
645 
646 extern N_gradient_2d *N_alloc_gradient_2d(void);
647 extern void N_free_gradient_2d(N_gradient_2d * grad);
648 extern N_gradient_2d *N_create_gradient_2d(double NC, double SC, double WC,
649  double EC);
650 extern int N_copy_gradient_2d(N_gradient_2d * source, N_gradient_2d * target);
652  N_gradient_2d * gradient, int col,
653  int row);
654 
655 extern N_gradient_3d *N_alloc_gradient_3d(void);
656 extern void N_free_gradient_3d(N_gradient_3d * grad);
657 extern N_gradient_3d *N_create_gradient_3d(double NC, double SC, double WC,
658  double EC, double TC, double BC);
659 extern int N_copy_gradient_3d(N_gradient_3d * source, N_gradient_3d * target);
661  N_gradient_3d * gradient, int col,
662  int row, int depth);
663 
667  double NEN,
668  double WC,
669  double EC,
670  double SWS,
671  double SES);
673  N_gradient_neighbours_x * target);
674 
678  double NEE,
679  double NC,
680  double SC,
681  double SWW,
682  double SEE);
684  N_gradient_neighbours_y * target);
685 
689  double NZ,
690  double NEZ,
691  double WZ,
692  double CZ,
693  double EZ,
694  double SWZ,
695  double SZ,
696  double SEZ);
698  N_gradient_neighbours_z * target);
699 
706  N_gradient_neighbours_2d * target);
709  N_gradient_neighbours_2d * gradient,
710  int col, int row);
711 
712 
725  N_gradient_neighbours_3d * target);
726 
729 
730 
731 extern N_gradient_field_2d *N_alloc_gradient_field_2d(int cols, int rows);
732 extern void N_free_gradient_field_2d(N_gradient_field_2d * field);
734  N_gradient_field_2d * target);
736  N_array_2d * weight_x,
737  N_array_2d * weight_y,
738  N_geom_data * geom,
740  gradfield);
742  field, N_array_2d * x_comp,
743  N_array_2d * y_comp);
744 
747 
749  int depths);
750 extern void N_free_gradient_field_3d(N_gradient_field_3d * field);
752  N_gradient_field_3d * target);
754  N_array_3d * weight_x,
755  N_array_3d * weight_y,
756  N_array_3d * weight_z,
757  N_geom_data * geom,
759  gradfield);
761  field, N_array_3d * x_comp,
762  N_array_3d * y_comp,
763  N_array_3d * z_comp);
764 
765 #endif
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), a weighting factor N_array_2d named weight and the distance between two cells saved in the N_geom_data struct.
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_gradient_neighbours_x * xt
Definition: N_pde.h:606
N_array_2d * N_math_array_2d(N_array_2d *array1, N_array_2d *array2, N_array_2d *result, int type)
Performe calculations with two input arrays, the result is written to a third array.
N_gradient_neighbours_z * zt
Definition: N_pde.h:614
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.
void N_free_gradient_neighbours_y(N_gradient_neighbours_y *grad)
Free's a N_gradient_neighbours_y structure.
Definition: N_gradient.c:403
Matrix entries for a mass balance 5/7/9 star system.
Definition: N_pde.h:341
Gradient between the cell neighbours in X and Y direction.
Definition: N_pde.h:593
N_array_3d * N_math_array_3d(N_array_3d *array1, N_array_3d *array2, N_array_3d *result, int type)
Performe calculations with two input arrays, the result is written to a third array.
float b
Definition: named_colr.c:8
int N_solver_pcg(N_les *les, int maxit, double error, int prec)
The iterative preconditioned conjugate gradients solver for symmetric positive definite matrices...
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)
int rows
Definition: N_pde.h:181
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
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:140
Gradient between the cell neighbours in Z direction.
Definition: N_pde.h:585
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 quad
Definition: N_pde.h:104
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:311
N_les_callback_2d * N_alloc_les_callback_2d(void)
Allocate the structure holding the callback function.
double N_calc_geom_mean(double a, double b)
Calculate the geometrical mean of values a and b.
Definition: N_tools.c:76
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:1018
void N_free_spvector(N_spvector *vector)
Release the memory of the sparse vector.
Definition: N_les.c:363
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:408
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:115
void N_print_les(N_les *les)
prints the linear equation system to stdout
Definition: N_les.c:312
int dim
Definition: N_pde.h:134
string name
Definition: render.py:1304
void N_free_gradient_neighbours_x(N_gradient_neighbours_x *grad)
Free's a N_gradient_neighbours_x structure.
Definition: N_gradient.c:307
#define min(x, y)
Definition: draw2.c:68
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:71
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:522
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
Gradient between the cells in X and Y direction.
Definition: N_pde.h:509
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:1018
N_data_star * N_callback_template_2d(void *data, N_geom_data *geom, int col, int row)
A callback template creates a 9 point star structure.
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:181
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:178
N_gradient_neighbours_y * y
Definition: N_pde.h:597
N_data_star * N_alloc_7star(void)
allocate a 7 point star data structure
void N_print_gradient_field_3d_info(N_gradient_field_3d *field)
Print gradient field information to stdout.
Definition: N_gradient.c:1090
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:96
int N_get_array_3d_type(N_array_3d *array3d)
Return the data type of the N_array_3d.
Definition: N_arrays.c:806
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.
N_STD_OPT
Standard options of the gpde library.
Definition: N_pde.h:445
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.
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:343
N_gradient_neighbours_z * N_alloc_gradient_neighbours_z(void)
Allocate a N_gradient_neighbours_z structure.
Definition: N_gradient.c:481
#define NE
Definition: dataquad.h:18
void N_free_array_2d(N_array_2d *data_array)
Release the memory of a N_array_2d structure.
Definition: N_arrays.c:127
int cols
Definition: N_pde.h:81
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:551
void N_calc_gradient_field_3d_stats(N_gradient_field_3d *field)
Calculate basic statistics of a gradient field.
N_gradient_neighbours_y * N_alloc_gradient_neighbours_y(void)
Allocate a N_gradient_neighbours_y structure.
Definition: N_gradient.c:386
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:826
void N_free_les(N_les *les)
Release the memory of the linear equation system.
Definition: N_les.c:387
double WC
Definition: N_pde.h:520
int N_copy_gradient_2d(N_gradient_2d *source, N_gradient_2d *target)
copy a N_gradient_2d structure
Definition: N_gradient.c:85
double N_calc_quad_mean(double a, double b)
Calculate the quadratic mean of values a and b.
Definition: N_tools.c:169
N_data_star * N_alloc_27star(void)
allocate a 27 point star data structure
void N_free_geom_data(N_geom_data *geodata)
Release memory of a pde geometry data structure.
Definition: N_geom.c:50
#define NW
Definition: dataquad.h:17
double * x
Definition: N_pde.h:98
int N_solver_cholesky(N_les *les)
The choleksy decomposition solver for quardatic, symmetric positiv definite matrices.
Definition: N_solvers.c:152
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_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:825
int y
Definition: plot.c:34
callback structure for 2d matrix assembling
Definition: N_pde.h:363
int depths
Definition: N_pde.h:142
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.
#define max(x, y)
Definition: draw2.c:69
N_array_3d * y_array
Definition: N_pde.h:637
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:329
N_geom_data * N_init_geom_data_3d(G3D_Region *region3d, N_geom_data *geodata)
Initiate a pde geometry data structure with a 3d region.
Definition: N_geom.c:73
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:198
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...
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:103
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 ...
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:1069
int offset
Definition: N_pde.h:225
int N_les_pivot_create(N_les *les)
Optimize the structure of the linear equation system with a common pivoting strategy.
Definition: N_les_pivot.c:47
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:408
int count
Definition: N_pde.h:344
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:1172
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...
double dy
Definition: N_pde.h:137
The row vector of the sparse matrix.
Definition: N_pde.h:79
N_gradient_neighbours_z * zb
Definition: N_pde.h:615
int offset
Definition: N_pde.h:183
N_les * N_alloc_les_Ax_b(int rows, int type)
Allocate memory for a quadratic linear equation system which includes the Matrix A, vector x and vector b.
Definition: N_les.c:165
tuple size
value.Bind(wx.EVT_TEXT, self.OnVolumeIsosurfMap)
Definition: tools.py:2334
double W_B
Definition: N_pde.h:349
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:248
CELL * cell_array
Definition: N_pde.h:184
#define SW
Definition: dataquad.h:19
void N_free_gradient_neighbours_z(N_gradient_neighbours_z *grad)
Free's a N_gradient_neighbours_z structure.
Definition: N_gradient.c:498
#define D
Definition: gis/intersect.c:74
#define C
Definition: intr_char.c:17
Geometric information about the structured grid.
Definition: N_pde.h:130
N_array_2d * y_array
Definition: N_pde.h:625
double dz
Definition: N_pde.h:138
int N_add_spvector_to_les(N_les *les, N_spvector *vector, int row)
Adds a sparse vector to a sparse linear equation system at position row.
Definition: N_les.c:263
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:46
tuple data
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:192
int rows
Definition: N_pde.h:143
N_data_star * N_alloc_5star(void)
allocate a 5 point star data structure
N_gradient_neighbours_x * xb
Definition: N_pde.h:608
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)
double ** A
Definition: N_pde.h:100
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:55
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:187
N_gradient_neighbours_y * yb
Definition: N_pde.h:612
int * index
Definition: N_pde.h:83
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:634
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 postion col, row is of type null, otherwise 0...
Definition: N_arrays.c:228
void N_print_array_3d(N_array_3d *data)
Write info and content of the array data to stdout.
Definition: N_arrays.c:1220
Gradient between the cell neighbours in X, Y and Z direction.
Definition: N_pde.h:603
void N_calc_gradient_field_2d_stats(N_gradient_field_2d *field)
Calculate basic statistics of a gradient field.
int N_get_array_2d_type(N_array_2d *array2d)
Return the data type of the N_array_2d struct.
Definition: N_arrays.c:161
N_array_3d * z_array
Definition: N_pde.h:638
Gradient between the cell neighbours in X direction.
Definition: N_pde.h:569
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:987
double W_T
Definition: N_pde.h:347
void N_free_gradient_neighbours_3d(N_gradient_neighbours_3d *grad)
Free's a N_gradient_neighbours_3d structure.
Definition: N_gradient.c:796
void N_sparse_matrix_vector_product(N_les *les, double *source, double *result)
Calculates the matrix - vector product of sparse matrix L->Asp and vector x.
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:189
N_spvector * N_alloc_spvector(int cols)
Allocate memory for a sparse vector.
Definition: N_les.c:29
double N_exp_upwinding(double sprod, double distance, double D)
exponential upwinding stabilization algorithm
Definition: N_upwind.c:63
N_gradient_neighbours_y * yt
Definition: N_pde.h:610
FCELL * fcell_array
Definition: N_pde.h:185
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:581
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:608
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)...
double Az
Definition: N_pde.h:140
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:702
double N_norm_array_3d(N_array_3d *array1, N_array_3d *array2, int type)
Calculate the norm of the two input arrays.
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...
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
void N_free_gradient_field_2d(N_gradient_field_2d *field)
Free's a N_gradient_neighbours_2d structure.
Definition: N_gradient.c:944
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:360
char * value
Definition: env.c:30
DCELL * dcell_array
Definition: N_pde.h:186
void N_free_array_3d(N_array_3d *data_array)
Release the memory of a N_array_3d.
Definition: N_arrays.c:777
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)
int N_solver_jacobi(N_les *L, int maxit, double sor, double error)
The iterative jacobian solver for regular matrices.
int N_solver_gauss(N_les *les)
The gauss elimination solver for quardatic matrices.
Definition: N_solvers.c:46
int cols
Definition: N_pde.h:144
Gradient between the cells in X, Y and Z direction.
Definition: N_pde.h:517
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:958
int type
Definition: N_pde.h:343
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:663
#define SE
Definition: dataquad.h:20
double N_calc_arith_mean(double a, double b)
Calculate the arithmetic mean of values a and b.
Definition: N_tools.c:33
float * fcell_array
Definition: N_pde.h:226
N_gradient_3d * N_alloc_gradient_3d(void)
Allocate a N_gradient_3d structure.
Definition: N_gradient.c:151
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:521
N_data_star * N_callback_template_3d(void *data, N_geom_data *geom, int col, int row, int depth)
A callback template creates a 7 point star structure.
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:557
int N_solver_cg(N_les *les, int maxit, double error)
The iterative conjugate gradients solver for symmetric positive definite matrices.
int planimetric
Definition: N_pde.h:132
N_les_callback_3d * N_alloc_les_callback_3d(void)
Allocate the structure holding the callback function.
double * dcell_array
Definition: N_pde.h:227
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:878
double * values
Definition: N_pde.h:82
void N_free_gradient_neighbours_2d(N_gradient_neighbours_2d *grad)
Free's a N_gradient_neighbours_2d structure.
Definition: N_gradient.c:608
int type
Definition: N_pde.h:180
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:455
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:723
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:630
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:1145
int cols
Definition: N_pde.h:103
tuple cols
N_gradient_neighbours_x * N_alloc_gradient_neighbours_x(void)
Allocate a N_gradient_neighbours_x structure.
Definition: N_gradient.c:290
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:455
void N_set_les_callback_3d_func(N_les_callback_3d *data, N_data_star *(*callback_func_3d)())
Set the callback function which is called while assembling the les in 3d.
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:135
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:1073
double * area
Definition: N_pde.h:133
N_array_3d * x_array
Definition: N_pde.h:636
double WC
Definition: N_pde.h:512
int N_convert_array_3d_null_to_zero(N_array_3d *a)
Convert all null values to zero values.
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:875
void N_free_gradient_field_3d(N_gradient_field_3d *field)
Free's a N_gradient_neighbours_3d structure.
Definition: N_gradient.c:1045
N_gradient_neighbours_x * xc
Definition: N_pde.h:607
int N_solver_lu(N_les *les)
The LU solver for quardatic matrices.
Definition: N_solvers.c:82
#define N
Definition: inverse.c:8
N_gradient_field_2d * N_alloc_gradient_field_2d(int cols, int rows)
Allocate a N_gradient_field_2d.
Definition: N_gradient.c:920
void N_free_gradient_2d(N_gradient_2d *grad)
Free's a N_gradient_2d structure.
Definition: N_gradient.c:42
void N_print_array_3d_info(N_array_3d *data)
Write the info of the array to stdout.
Definition: N_arrays.c:1194
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:424
double * b
Definition: N_pde.h:99
N_gradient_neighbours_y * yc
Definition: N_pde.h:611
double N_calc_harmonic_mean(double a, double b)
Calculate the harmonical mean of values a and b.
Definition: N_tools.c:119
int N_solver_bicgstab(N_les *les, int maxit, double error)
The iterative biconjugate gradients solver with stabilization for unsymmetric non-definite matrices...
int N_copy_gradient_3d(N_gradient_3d *source, N_gradient_3d *target)
copy a N_gradient_3d structure
Definition: N_gradient.c:214
int N_convert_array_2d_null_to_zero(N_array_2d *a)
Convert all null values to zero values.
callback structure for 3d matrix assembling
Definition: N_pde.h:355
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:87
Gradient between the cell neighbours in Y direction.
Definition: N_pde.h:577
double W
Definition: N_pde.h:345
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:276
double dx
Definition: N_pde.h:136
int rows
Definition: N_pde.h:102
void N_free_gradient_3d(N_gradient_3d *grad)
Free's a N_gradient_3d structure.
Definition: N_gradient.c:166
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:114
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
N_gradient_neighbours_3d * N_alloc_gradient_neighbours_3d(void)
Allocate a N_gradient_neighbours_3d structure.
Definition: N_gradient.c:769
N_gradient_neighbours_x * x
Definition: N_pde.h:596
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:966
N_data_star * N_alloc_9star(void)
allocate a 9 point star data structure
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...
int type
Definition: N_pde.h:222
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_alloc_les_A(int rows, int type)
Allocate memory for a quadratic linear equation system which includes the Matrix A.
Definition: N_les.c:150
int N_solver_SOR(N_les *L, int maxit, double sor, double error)
The iterative overrelaxed gauss seidel solver for regular matrices.
void N_print_gradient_field_2d_info(N_gradient_field_2d *field)
Print gradient field information to stdout.
Definition: N_gradient.c:986
int rows_intern
Definition: N_pde.h:224
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:375
void N_set_les_callback_2d_func(N_les_callback_2d *data, N_data_star *(*callback_func_2d)())
Set the callback function which is called while assembling the les in 2d.
SC
Definition: wxgui.py:32
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.
The linear equation system (les) structure.
Definition: N_pde.h:96
N_array_2d * x_array
Definition: N_pde.h:624
N_gradient_2d * N_alloc_gradient_2d(void)
Allocate a N_gradient_2d structure.
Definition: N_gradient.c:27
struct Option * N_define_standard_option(int opt)
Create standardised Option structure related to the gpde library.
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:30
int rows
Definition: N_pde.h:223
N_spvector ** Asp
Definition: N_pde.h:101
N_les * N_create_diag_precond_matrix(N_les *les, int prec)
Compute a diagonal preconditioning matrix for krylov space solver.
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:60
int type
Definition: N_pde.h:105
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:69
double N_full_upwinding(double sprod, double distance, double D)
full upwinding stabilization algorithm
Definition: N_upwind.c:33
int rows_intern
Definition: N_pde.h:182
N_les * N_alloc_les(int rows, int type)
Allocate memory for a quadratic linear equation system which includes the Matrix A, vector x and vector b.
Definition: N_les.c:120
void N_matrix_vector_product(N_les *les, double *source, double *result)
Calculates the matrix - vector product of matrix L->A and vector x.
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:53
N_gradient_neighbours_2d * N_alloc_gradient_neighbours_2d(void)
Allocate a N_gradient_neighbours_2d structure.
Definition: N_gradient.c:587