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/raster3d.h> 20 #include <grass/glocale.h> 21 #include <grass/gmath.h> 22 23 #ifndef _N_PDE_H_ 24 #define _N_PDE_H_ 25 26 #define N_NORMAL_LES 0 27 #define N_SPARSE_LES 1 28 /*! 29 * Boundary conditions for cells 30 */ 31 #define N_CELL_INACTIVE 0 32 #define N_CELL_ACTIVE 1 33 #define N_CELL_DIRICHLET 2 34 #define N_CELL_TRANSMISSION 3 35 /*! 36 * \brief the maximum number of available cell states (eg: boundary condition, 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 3 /* array division, if div with 0 the NULL value is set */ 52 53 #define N_UPWIND_FULL 0 /*full upwinding stabilization */ 54 #define N_UPWIND_EXP 1 /*exponential upwinding stabilization */ 55 #define N_UPWIND_WEIGHT 2 /*weighted upwinding stabilization */ 56 57 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 { 74 double *x; /*the value vector */ 75 double *b; /*the right side of Ax = b */ 76 double **A; /*the normal quadratic matrix */ 77 G_math_spvector **Asp; /*the sparse matrix */ 78 int rows; /*number of rows */ 79 int cols; /*number of cols */ 80 int quad; /*is the matrix quadratic (1-quadratic, 0 not) */ 81 int type; /*the type of the les, normal == 0, sparse == 1 */ 82 } N_les; 83 84 extern N_les *N_alloc_les_param(int cols, int rows, int type, int param); 85 extern N_les *N_alloc_les(int rows, int type); 86 extern N_les *N_alloc_les_A(int rows, int type); 87 extern N_les *N_alloc_les_Ax(int rows, int type); 88 extern N_les *N_alloc_les_Ax_b(int rows, int type); 89 extern N_les *N_alloc_nquad_les(int cols, int rows, int type); 90 extern N_les *N_alloc_nquad_les_A(int cols, int rows, int type); 91 extern N_les *N_alloc_nquad_les_Ax(int cols, int rows, int type); 92 extern N_les *N_alloc_nquad_les_Ax_b(int cols, int rows, int type); 93 extern void N_print_les(N_les * les); 94 extern void N_free_les(N_les * les); 95 96 /* *************************************************************** */ 97 /* *************** GEOMETRY INFORMATION ************************** */ 98 /* *************************************************************** */ 99 100 /*! 101 * \brief Geometric information about the structured grid 102 * */ 103 typedef struct 104 { 105 int planimetric; /*If the projection is not planimetric (0), the array calculation is different for each row */ 106 double *area; /* the vector of area values for non-planimetric projection 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); 123 extern N_geom_data *N_init_geom_data_3d(RASTER3D_Region * region3d, N_geom_data * geodata); 124 extern N_geom_data *N_init_geom_data_2d(struct Cell_head *region, N_geom_data * geodata); 125 extern double N_get_geom_data_area_of_cell(N_geom_data * geom, int row); 126 127 /* *************************************************************** */ 128 /* *************** READING RASTER AND VOLUME DATA **************** */ 129 /* *************************************************************** */ 130 131 typedef struct 132 { 133 int type; /* which raster type CELL_TYPE, FCELL_TYPE, DCELL_TYPE */ 134 int rows, cols; 135 int rows_intern, cols_intern; 136 int offset; /*number of cols/rows offset at each boundary */ 137 CELL *cell_array; /*The data is stored in an one dimensional array internally */ 138 FCELL *fcell_array; /*The data is stored in an one dimensional array internally */ 139 DCELL *dcell_array; /*The data is stored in an one dimensional array internally */ 140 } N_array_2d; 141 142 extern N_array_2d *N_alloc_array_2d(int cols, int rows, int offset, int type); 143 extern void N_free_array_2d(N_array_2d * data_array); 144 extern int N_get_array_2d_type(N_array_2d * array2d); 145 extern void N_get_array_2d_value(N_array_2d * array2d, int col, int row, void *value); 146 extern CELL N_get_array_2d_c_value(N_array_2d * array2d, int col, int row); 147 extern FCELL N_get_array_2d_f_value(N_array_2d * array2d, int col, int row); 148 extern DCELL N_get_array_2d_d_value(N_array_2d * array2d, int col, int row); 149 extern void N_put_array_2d_value(N_array_2d * array2d, int col, int row, char *value); 150 extern void N_put_array_2d_c_value(N_array_2d * array2d, int col, int row, CELL value); 151 extern void N_put_array_2d_f_value(N_array_2d * array2d, int col, int row, FCELL value); 152 extern void N_put_array_2d_d_value(N_array_2d * array2d, int col, int row, DCELL value); 153 extern int N_is_array_2d_value_null(N_array_2d * array2d, int col, int row); 154 extern void N_put_array_2d_value_null(N_array_2d * array2d, int col, int row); 155 extern void N_print_array_2d(N_array_2d * data); 156 extern void N_print_array_2d_info(N_array_2d * data); 157 extern void N_copy_array_2d(N_array_2d * source, N_array_2d * target); 158 extern double N_norm_array_2d(N_array_2d * array1, N_array_2d * array2, int type); 159 extern N_array_2d *N_math_array_2d(N_array_2d * array1, N_array_2d * array2, N_array_2d * result, int type); 160 extern int N_convert_array_2d_null_to_zero(N_array_2d * a); 161 extern N_array_2d *N_read_rast_to_array_2d(char *name, N_array_2d * array); 162 extern void N_write_array_2d_to_rast(N_array_2d * array, char *name); 163 extern void N_calc_array_2d_stats(N_array_2d * a, double *min, double *max, double *sum, int *nonzero, int withoffset); 164 165 typedef struct 166 { 167 int type; /* which raster type FCELL_TYPE, DCELL_TYPE */ 168 int rows, cols, depths; 169 int rows_intern, cols_intern, depths_intern; 170 int offset; /*number of cols/rows/depths offset at each boundary */ 171 float *fcell_array; /*The data is stored in an one dimensional array internally */ 172 double *dcell_array; /*The data is stored in an one dimensional array internally */ 173 } N_array_3d; 174 175 extern N_array_3d *N_alloc_array_3d(int cols, int rows, int depths, int offset, int type); 176 extern void N_free_array_3d(N_array_3d * data_array); 177 extern int N_get_array_3d_type(N_array_3d * array3d); 178 extern void N_get_array_3d_value(N_array_3d * array3d, int col, int row, int depth, void *value); 179 extern float N_get_array_3d_f_value(N_array_3d * array3d, int col, int row, int depth); 180 extern double N_get_array_3d_d_value(N_array_3d * array3d, int col, int row, int depth); 181 extern void N_put_array_3d_value(N_array_3d * array3d, int col, int row, int depth, char *value); 182 extern void N_put_array_3d_f_value(N_array_3d * array3d, int col, int row, int depth, float value); 183 extern void N_put_array_3d_d_value(N_array_3d * array3d, int col, int row, int depth, double value); 184 extern int N_is_array_3d_value_null(N_array_3d * array3d, int col, int row, int depth); 185 extern void N_put_array_3d_value_null(N_array_3d * array3d, int col, int row, int depth); 186 extern void N_print_array_3d(N_array_3d * data); 187 extern void N_print_array_3d_info(N_array_3d * data); 188 extern void N_copy_array_3d(N_array_3d * source, N_array_3d * target); 189 extern double N_norm_array_3d(N_array_3d * array1, N_array_3d * array2, int type); 190 extern N_array_3d *N_math_array_3d(N_array_3d * array1, N_array_3d * array2, N_array_3d * result, int type); 191 extern int N_convert_array_3d_null_to_zero(N_array_3d * a); 192 extern N_array_3d *N_read_rast3d_to_array_3d(char *name, N_array_3d * array, int mask); 193 extern void N_write_array_3d_to_rast3d(N_array_3d * array, char *name, int mask); 194 extern void N_calc_array_3d_stats(N_array_3d * a, double *min, double *max, double *sum, int *nonzero, int withoffset); 195 196 /* *************************************************************** */ 197 /* *************** MATRIX ASSEMBLING METHODS ********************* */ 198 /* *************************************************************** */ 199 /*! 200 * \brief Matrix entries for a mass balance 5/7/9 star system 201 * 202 * Matrix entries for the mass balance of a 5 star system 203 * 204 * The entries are center, east, west, north, south and the 205 * right side vector b of Ax = b. This system is typically used in 2d. 206 207 \verbatim 208 N 209 | 210 W-- C --E 211 | 212 S 213 \endverbatim 214 215 * Matrix entries for the mass balance of a 7 star system 216 * 217 * The entries are center, east, west, north, south, top, bottom and the 218 * right side vector b of Ax = b. This system is typically used in 3d. 219 220 \verbatim 221 T N 222 |/ 223 W-- C --E 224 /| 225 S B 226 \endverbatim 227 228 * Matrix entries for the mass balance of a 9 star system 229 * 230 * The entries are center, east, west, north, south, north-east, south-east, 231 * north-wast, south-west and the 232 * right side vector b of Ax = b. This system is typically used in 2d. 233 234 \verbatim 235 NW N NE 236 \ | / 237 W-- C --E 238 / | \ 239 SW S SE 240 \endverbatim 241 242 * Matrix entries for the mass balance of a 27 star system 243 * 244 * The entries are center, east, west, north, south, north-east, south-east, 245 * north-wast, south-west, same for top and bottom and the 246 * right side vector b of Ax = b. This system is typically used in 2d. 247 248 \verbatim 249 top: 250 NW_T N_Z NE_T 251 \ | / 252 W_T-- T --E_T 253 / | \ 254 SW_T S_T SE_T 255 256 center: 257 NW N NE 258 \ | / 259 W-- C --E 260 / | \ 261 SW S SE 262 263 bottom: 264 NW_B N_B NE_B 265 \ | / 266 W_B-- B --E_B 267 / | \ 268 SW_B S_B SE_B 269 \endverbatim 270 271 */ 272 typedef struct 273 { 274 int type; 275 int count; 276 double C, W, E, N, S, NE, NW, SE, SW, V; 277 /*top part */ 278 double T, W_T, E_T, N_T, S_T, NE_T, NW_T, SE_T, SW_T; 279 /*bottom part */ 280 double B, W_B, E_B, N_B, S_B, NE_B, NW_B, SE_B, SW_B; 281 } N_data_star; 282 283 /*! 284 * \brief callback structure for 3d matrix assembling 285 * */ 286 typedef struct 287 { 288 N_data_star *(*callback) (); 289 } N_les_callback_3d; 290 291 /*! 292 * \brief callback structure for 2d matrix assembling 293 * */ 294 typedef struct 295 { 296 N_data_star *(*callback) (); 297 } N_les_callback_2d; 298 299 300 extern void N_set_les_callback_3d_func(N_les_callback_3d * data, N_data_star * (*callback_func_3d) ()); 301 extern void N_set_les_callback_2d_func(N_les_callback_2d * data, N_data_star * (*callback_func_2d) ()); 302 extern N_les_callback_3d *N_alloc_les_callback_3d(void); 303 extern N_les_callback_2d *N_alloc_les_callback_2d(void); 304 extern N_data_star *N_alloc_5star(void); 305 extern N_data_star *N_alloc_7star(void); 306 extern N_data_star *N_alloc_9star(void); 307 extern N_data_star *N_alloc_27star(void); 308 extern N_data_star *N_create_5star(double C, double W, double E, double N, 309 double S, double V); 310 extern N_data_star *N_create_7star(double C, double W, double E, double N, 311 double S, double T, double B, double V); 312 extern N_data_star *N_create_9star(double C, double W, double E, double N, 313 double S, double NW, double SW, double NE, 314 double SE, double V); 315 extern N_data_star *N_create_27star(double C, double W, double E, double N, 316 double S, double NW, double SW, double NE, 317 double SE, double T, double W_T, 318 double E_T, double N_T, double S_T, 319 double NW_T, double SW_T, double NE_T, 320 double SE_T, double B, double W_B, 321 double E_B, double N_B, double S_B, 322 double NW_B, double SW_B, double NE_B, 323 double SE_B, double V); 324 extern N_data_star *N_callback_template_3d(void *data, N_geom_data * geom, int col, int row, int depth); 325 extern N_data_star *N_callback_template_2d(void *data, N_geom_data * geom, int col, int row); 326 extern 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); 327 extern 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); 328 extern 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); 329 extern 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); 330 extern 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); 331 extern 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); 332 extern 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); 333 extern 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); 334 extern int N_les_pivot_create(N_les * les); 335 int N_les_integrate_dirichlet_2d(N_les * les, N_geom_data * geom, N_array_2d * status, N_array_2d * start_val); 336 int N_les_integrate_dirichlet_3d(N_les * les, N_geom_data * geom, N_array_3d * status, N_array_3d * start_val); 337 338 /* *************************************************************** */ 339 /* *************** GPDE STANDARD OPTIONS ************************* */ 340 /* *************************************************************** */ 341 342 /*! \brief Standard options of the gpde library 343 * */ 344 typedef enum 345 { 346 N_OPT_SOLVER_SYMM, /*! solver for symmetric, positive definite linear equation systems */ 347 N_OPT_SOLVER_UNSYMM, /*! solver for unsymmetric linear equation systems */ 348 N_OPT_MAX_ITERATIONS, /*! Maximum number of iteration used to solver the linear equation system */ 349 N_OPT_ITERATION_ERROR, /*! Error break criteria for the iterative solver (jacobi, sor, cg or bicgstab) */ 350 N_OPT_SOR_VALUE, /*! The relaxation parameter used by the jacobi and sor solver for speedup or stabilizing */ 351 N_OPT_CALC_TIME /*! The calculation time in seconds */ 352 } N_STD_OPT; 353 354 extern struct Option *N_define_standard_option(int opt); 355 356 /* *************************************************************** */ 357 /* *************** GPDE MATHEMATICAL TOOLS *********************** */ 358 /* *************************************************************** */ 359 360 extern double N_calc_arith_mean(double a, double b); 361 extern double N_calc_arith_mean_n(double *a, int size); 362 extern double N_calc_geom_mean(double a, double b); 363 extern double N_calc_geom_mean_n(double *a, int size); 364 extern double N_calc_harmonic_mean(double a, double b); 365 extern double N_calc_harmonic_mean_n(double *a, int size); 366 extern double N_calc_quad_mean(double a, double b); 367 extern double N_calc_quad_mean_n(double *a, int size); 368 369 /* *************************************************************** */ 370 /* *************** UPWIND STABILIZATION ALGORITHMS *************** */ 371 /* *************************************************************** */ 372 373 extern double N_full_upwinding(double sprod, double distance, double D); 374 extern double N_exp_upwinding(double sprod, double distance, double D); 375 376 377 /* *************************************************************** */ 378 /* *************** METHODS FOR GRADIENT CALCULATION ************** */ 379 /* *************************************************************** */ 380 /*! 381 \verbatim 382 383 ______________ 384 | | | | 385 | | | | 386 |----|-NC-|----| 387 | | | | 388 | WC EC | 389 | | | | 390 |----|-SC-|----| 391 | | | | 392 |____|____|____| 393 394 395 | / 396 TC NC 397 |/ 398 --WC-----EC-- 399 /| 400 SC BC 401 / | 402 403 \endverbatim 404 405 */ 406 407 /*! \brief Gradient between the cells in X and Y direction */ 408 typedef struct 409 { 410 411 double NC, SC, WC, EC; 412 413 } N_gradient_2d; 414 415 /*! \brief Gradient between the cells in X, Y and Z direction */ 416 typedef struct 417 { 418 419 double NC, SC, WC, EC, TC, BC; 420 421 } N_gradient_3d; 422 423 424 /*! 425 \verbatim 426 427 Gradient in X direction between the cell neighbours 428 ____ ____ ____ 429 | | | | 430 | NWN NEN | 431 |____|____|____| 432 | | | | 433 | WN EN | 434 |____|____|____| 435 | | | | 436 | SWS SES | 437 |____|____|____| 438 439 Gradient in Y direction between the cell neighbours 440 ______________ 441 | | | | 442 | | | | 443 |NWW-|-NC-|-NEE| 444 | | | | 445 | | | | 446 |SWW-|-SC-|-SEE| 447 | | | | 448 |____|____|____| 449 450 Gradient in Z direction between the cell neighbours 451 /______________/ 452 /| | | | 453 | NWZ| NZ | NEZ| 454 |____|____|____| 455 /| | | | 456 | WZ | CZ | EZ | 457 |____|____|____| 458 /| | | | 459 | SWZ| SZ | SEZ| 460 |____|____|____| 461 /____/____/____/ 462 463 464 \endverbatim 465 */ 466 467 /*! \brief Gradient between the cell neighbours in X direction */ 468 typedef struct 469 { 470 471 double NWN, NEN, WC, EC, SWS, SES; 472 473 } N_gradient_neighbours_x; 474 475 /*! \brief Gradient between the cell neighbours in Y direction */ 476 typedef struct 477 { 478 479 double NWW, NEE, NC, SC, SWW, SEE; 480 481 } N_gradient_neighbours_y; 482 483 /*! \brief Gradient between the cell neighbours in Z direction */ 484 typedef struct 485 { 486 487 double NWZ, NZ, NEZ, WZ, CZ, EZ, SWZ, SZ, SEZ; 488 489 } N_gradient_neighbours_z; 490 491 /*! \brief Gradient between the cell neighbours in X and Y direction */ 492 typedef struct 493 { 494 495 N_gradient_neighbours_x *x; 496 N_gradient_neighbours_y *y; 497 498 } N_gradient_neighbours_2d; 499 500 501 /*! \brief Gradient between the cell neighbours in X, Y and Z direction */ 502 typedef struct 503 { 504 505 N_gradient_neighbours_x *xt; /*top values */ 506 N_gradient_neighbours_x *xc; /*center values */ 507 N_gradient_neighbours_x *xb; /*bottom values */ 508 509 N_gradient_neighbours_y *yt; /*top values */ 510 N_gradient_neighbours_y *yc; /*center values */ 511 N_gradient_neighbours_y *yb; /*bottom values */ 512 513 N_gradient_neighbours_z *zt; /*top-center values */ 514 N_gradient_neighbours_z *zb; /*bottom-center values */ 515 516 } N_gradient_neighbours_3d; 517 518 519 /*! Two dimensional gradient field */ 520 typedef struct 521 { 522 523 N_array_2d *x_array; 524 N_array_2d *y_array; 525 int cols, rows; 526 double min, max, mean, sum; 527 int nonull; 528 529 } N_gradient_field_2d; 530 531 /*! Three dimensional gradient field */ 532 typedef struct 533 { 534 535 N_array_3d *x_array; 536 N_array_3d *y_array; 537 N_array_3d *z_array; 538 int cols, rows, depths; 539 double min, max, mean, sum; 540 int nonull; 541 542 } N_gradient_field_3d; 543 544 545 extern N_gradient_2d *N_alloc_gradient_2d(void); 546 extern void N_free_gradient_2d(N_gradient_2d * grad); 547 extern N_gradient_2d *N_create_gradient_2d(double NC, double SC, double WC, double EC); 548 extern int N_copy_gradient_2d(N_gradient_2d * source, N_gradient_2d * target); 549 extern N_gradient_2d *N_get_gradient_2d(N_gradient_field_2d * field, N_gradient_2d * gradient, int col, int row); 550 extern N_gradient_3d *N_alloc_gradient_3d(void); 551 extern void N_free_gradient_3d(N_gradient_3d * grad); 552 extern N_gradient_3d *N_create_gradient_3d(double NC, double SC, double WC, double EC, double TC, double BC); 553 extern int N_copy_gradient_3d(N_gradient_3d * source, N_gradient_3d * target); 554 extern N_gradient_3d *N_get_gradient_3d(N_gradient_field_3d * field, N_gradient_3d * gradient, int col, int row, int depth); 555 extern N_gradient_neighbours_x *N_alloc_gradient_neighbours_x(void); 556 extern void N_free_gradient_neighbours_x(N_gradient_neighbours_x * grad); 557 extern N_gradient_neighbours_x *N_create_gradient_neighbours_x(double NWN, 558 double NEN, 559 double WC, 560 double EC, 561 double SWS, 562 double SES); 563 extern int N_copy_gradient_neighbours_x(N_gradient_neighbours_x * source, N_gradient_neighbours_x * target); 564 extern N_gradient_neighbours_y *N_alloc_gradient_neighbours_y(void); 565 extern void N_free_gradient_neighbours_y(N_gradient_neighbours_y * grad); 566 extern N_gradient_neighbours_y *N_create_gradient_neighbours_y(double NWW, 567 double NEE, 568 double NC, 569 double SC, 570 double SWW, 571 double SEE); 572 extern int N_copy_gradient_neighbours_y(N_gradient_neighbours_y * source, N_gradient_neighbours_y * target); 573 extern N_gradient_neighbours_z *N_alloc_gradient_neighbours_z(void); 574 extern void N_free_gradient_neighbours_z(N_gradient_neighbours_z * grad); 575 extern N_gradient_neighbours_z *N_create_gradient_neighbours_z(double NWZ, 576 double NZ, 577 double NEZ, 578 double WZ, 579 double CZ, 580 double EZ, 581 double SWZ, 582 double SZ, 583 double SEZ); 584 extern int N_copy_gradient_neighbours_z(N_gradient_neighbours_z * source, N_gradient_neighbours_z * target); 585 extern N_gradient_neighbours_2d *N_alloc_gradient_neighbours_2d(void); 586 extern void N_free_gradient_neighbours_2d(N_gradient_neighbours_2d * grad); 587 extern N_gradient_neighbours_2d * N_create_gradient_neighbours_2d(N_gradient_neighbours_x * x, N_gradient_neighbours_y * y); 588 extern int N_copy_gradient_neighbours_2d(N_gradient_neighbours_2d * source, N_gradient_neighbours_2d * target); 589 extern N_gradient_neighbours_2d * N_get_gradient_neighbours_2d(N_gradient_field_2d * field, N_gradient_neighbours_2d * gradient, int col, int row); 590 extern N_gradient_neighbours_3d *N_alloc_gradient_neighbours_3d(void); 591 extern void N_free_gradient_neighbours_3d(N_gradient_neighbours_3d * grad); 592 extern N_gradient_neighbours_3d 593 * N_create_gradient_neighbours_3d(N_gradient_neighbours_x * xt, 594 N_gradient_neighbours_x * xc, 595 N_gradient_neighbours_x * xb, 596 N_gradient_neighbours_y * yt, 597 N_gradient_neighbours_y * yc, 598 N_gradient_neighbours_y * yb, 599 N_gradient_neighbours_z * zt, 600 N_gradient_neighbours_z * zb); 601 extern int N_copy_gradient_neighbours_3d(N_gradient_neighbours_3d * source, N_gradient_neighbours_3d * target); 602 extern void N_print_gradient_field_2d_info(N_gradient_field_2d * field); 603 extern void N_calc_gradient_field_2d_stats(N_gradient_field_2d * field); 604 extern N_gradient_field_2d *N_alloc_gradient_field_2d(int cols, int rows); 605 extern void N_free_gradient_field_2d(N_gradient_field_2d * field); 606 extern int N_copy_gradient_field_2d(N_gradient_field_2d * source, N_gradient_field_2d * target); 607 extern N_gradient_field_2d *N_compute_gradient_field_2d(N_array_2d * pot, 608 N_array_2d * weight_x, 609 N_array_2d * weight_y, 610 N_geom_data * geom, 611 N_gradient_field_2d * 612 gradfield); 613 extern void N_compute_gradient_field_components_2d(N_gradient_field_2d * field, N_array_2d * x_comp, N_array_2d * y_comp); 614 extern void N_print_gradient_field_3d_info(N_gradient_field_3d * field); 615 extern void N_calc_gradient_field_3d_stats(N_gradient_field_3d * field); 616 extern N_gradient_field_3d *N_alloc_gradient_field_3d(int cols, int rows, int depths); 617 extern void N_free_gradient_field_3d(N_gradient_field_3d * field); 618 extern int N_copy_gradient_field_3d(N_gradient_field_3d * source, N_gradient_field_3d * target); 619 extern N_gradient_field_3d *N_compute_gradient_field_3d(N_array_3d * pot, 620 N_array_3d * weight_x, 621 N_array_3d * weight_y, 622 N_array_3d * weight_z, 623 N_geom_data * geom, 624 N_gradient_field_3d * 625 gradfield); 626 extern 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); 627 628 #endif 629