GRASS Programmer's Manual
6.4.2(2012)
|
00001 00002 /***************************************************************************** 00003 * 00004 * MODULE: Grass PDE Numerical Library 00005 * AUTHOR(S): Soeren Gebbert, Berlin (GER) Dec 2006 00006 * soerengebbert <at> gmx <dot> de 00007 * 00008 * PURPOSE: This file contains definitions of variables and data types 00009 * 00010 * COPYRIGHT: (C) 2000 by the GRASS Development Team 00011 * 00012 * This program is free software under the GNU General Public 00013 * License (>=v2). Read the file COPYING that comes with GRASS 00014 * for details. 00015 * 00016 *****************************************************************************/ 00017 00018 #include <grass/gis.h> 00019 #include <grass/G3d.h> 00020 #include <grass/glocale.h> 00021 00022 #ifndef _N_PDE_H_ 00023 #define _N_PDE_H_ 00024 00025 /*solver names */ 00026 #define N_SOLVER_DIRECT_GAUSS "gauss" 00027 #define N_SOLVER_DIRECT_LU "lu" 00028 #define N_SOLVER_DIRECT_CHOLESKY "cholesky" 00029 #define N_SOLVER_ITERATIVE_JACOBI "jacobi" 00030 #define N_SOLVER_ITERATIVE_SOR "sor" 00031 #define N_SOLVER_ITERATIVE_CG "cg" 00032 #define N_SOLVER_ITERATIVE_PCG "pcg" 00033 #define N_SOLVER_ITERATIVE_BICGSTAB "bicgstab" 00034 00035 /*preconditioner */ 00036 #define N_DIAGONAL_PRECONDITION 1 00037 #define N_ROWSCALE_ABSSUMNORM_PRECONDITION 2 00038 #define N_ROWSCALE_EUKLIDNORM_PRECONDITION 3 00039 #define N_ROWSCALE_MAXNORM_PRECONDITION 4 00040 00041 #define N_NORMAL_LES 0 00042 #define N_SPARSE_LES 1 00043 00044 #define N_CELL_INACTIVE 0 00045 #define N_CELL_ACTIVE 1 00046 #define N_CELL_DIRICHLET 2 00047 #define N_CELL_TRANSMISSION 3 00048 00051 #define N_MAX_CELL_STATE 20 00052 00053 #define N_5_POINT_STAR 0 00054 #define N_7_POINT_STAR 1 00055 #define N_9_POINT_STAR 2 00056 #define N_27_POINT_STAR 3 00057 00058 #define N_MAXIMUM_NORM 0 00059 #define N_EUKLID_NORM 1 00060 00061 #define N_ARRAY_SUM 0 /* summ two arrays */ 00062 #define N_ARRAY_DIF 1 /* calc the difference between two arrays */ 00063 #define N_ARRAY_MUL 2 /* multiply two arrays */ 00064 #define N_ARRAY_DIV 3 /* array division, if div with 0 the NULL value is set */ 00065 00066 #define N_UPWIND_FULL 0 /*full upwinding stabilization */ 00067 #define N_UPWIND_EXP 1 /*exponential upwinding stabilization */ 00068 #define N_UPWIND_WEIGHT 2 /*weighted upwinding stabilization */ 00069 00070 00071 00072 /* *************************************************************** */ 00073 /* *************** LINEARE EQUATION SYSTEM PART ****************** */ 00074 /* *************************************************************** */ 00075 00079 typedef struct 00080 { 00081 int cols; /*Number of entries */ 00082 double *values; /*The non null values of the row */ 00083 int *index; /*the index number */ 00084 } N_spvector; 00085 00086 00096 typedef struct 00097 { 00098 double *x; /*the value vector */ 00099 double *b; /*the right side of Ax = b */ 00100 double **A; /*the normal quadratic matrix */ 00101 N_spvector **Asp; /*the sparse matrix */ 00102 int rows; /*number of rows */ 00103 int cols; /*number of cols */ 00104 int quad; /*is the matrix quadratic (1-quadratic, 0 not) */ 00105 int type; /*the type of the les, normal == 0, sparse == 1 */ 00106 } N_les; 00107 00108 extern N_spvector *N_alloc_spvector(int cols); 00109 extern N_les *N_alloc_les_param(int cols, int rows, int type, int param); 00110 extern N_les *N_alloc_les(int rows, int type); 00111 extern N_les *N_alloc_les_A(int rows, int type); 00112 extern N_les *N_alloc_les_Ax(int rows, int type); 00113 extern N_les *N_alloc_les_Ax_b(int rows, int type); 00114 extern N_les *N_alloc_nquad_les(int cols, int rows, int type); 00115 extern N_les *N_alloc_nquad_les_A(int cols, int rows, int type); 00116 extern N_les *N_alloc_nquad_les_Ax(int cols, int rows, int type); 00117 extern N_les *N_alloc_nquad_les_Ax_b(int cols, int rows, int type); 00118 extern void N_print_les(N_les * les); 00119 extern int N_add_spvector_to_les(N_les * les, N_spvector * vector, int row); 00120 extern void N_free_spvector(N_spvector * vector); 00121 extern void N_free_les(N_les * les); 00122 00123 /* *************************************************************** */ 00124 /* *************** GEOMETRY INFORMATION ************************** */ 00125 /* *************************************************************** */ 00126 00130 typedef struct 00131 { 00132 int planimetric; /*If the projection is not planimetric (0), the array calculation is different for each row */ 00133 double *area; /* the vector of area values for non-planimetric projection for each row */ 00134 int dim; /* 2 or 3 */ 00135 00136 double dx; 00137 double dy; 00138 double dz; 00139 00140 double Az; 00141 00142 int depths; 00143 int rows; 00144 int cols; 00145 00146 } N_geom_data; 00147 00148 extern N_geom_data *N_alloc_geom_data(void); 00149 extern void N_free_geom_data(N_geom_data * geodata); 00150 extern N_geom_data *N_init_geom_data_3d(G3D_Region * region3d, 00151 N_geom_data * geodata); 00152 extern N_geom_data *N_init_geom_data_2d(struct Cell_head *region, 00153 N_geom_data * geodata); 00154 extern double N_get_geom_data_area_of_cell(N_geom_data * geom, int row); 00155 00156 00157 /* *************************************************************** */ 00158 /* *************** LINEARE EQUATION SOLVER PART ****************** */ 00159 /* *************************************************************** */ 00160 extern int N_solver_gauss(N_les * les); 00161 extern int N_solver_lu(N_les * les); 00162 extern int N_solver_cholesky(N_les * les); 00163 extern int N_solver_jacobi(N_les * L, int maxit, double sor, double error); 00164 extern int N_solver_SOR(N_les * L, int maxit, double sor, double error); 00165 extern int N_solver_cg(N_les * les, int maxit, double error); 00166 extern int N_solver_pcg(N_les * les, int maxit, double error, int prec); 00167 extern int N_solver_bicgstab(N_les * les, int maxit, double error); 00168 extern void N_matrix_vector_product(N_les * les, double *source, 00169 double *result); 00170 extern void N_sparse_matrix_vector_product(N_les * les, double *source, 00171 double *result); 00172 extern N_les *N_create_diag_precond_matrix(N_les * les, int prec); 00173 00174 /* *************************************************************** */ 00175 /* *************** READING RASTER AND VOLUME DATA **************** */ 00176 /* *************************************************************** */ 00177 00178 typedef struct 00179 { 00180 int type; /* which raster type CELL_TYPE, FCELL_TYPE, DCELL_TYPE */ 00181 int rows, cols; 00182 int rows_intern, cols_intern; 00183 int offset; /*number of cols/rows offset at each boundary */ 00184 CELL *cell_array; /*The data is stored in an one dimensional array internally */ 00185 FCELL *fcell_array; /*The data is stored in an one dimensional array internally */ 00186 DCELL *dcell_array; /*The data is stored in an one dimensional array internally */ 00187 } N_array_2d; 00188 00189 extern N_array_2d *N_alloc_array_2d(int cols, int rows, int offset, int type); 00190 extern void N_free_array_2d(N_array_2d * data_array); 00191 extern int N_get_array_2d_type(N_array_2d * array2d); 00192 extern void N_get_array_2d_value(N_array_2d * array2d, int col, int row, 00193 void *value); 00194 extern CELL N_get_array_2d_c_value(N_array_2d * array2d, int col, int row); 00195 extern FCELL N_get_array_2d_f_value(N_array_2d * array2d, int col, int row); 00196 extern DCELL N_get_array_2d_d_value(N_array_2d * array2d, int col, int row); 00197 extern void N_put_array_2d_value(N_array_2d * array2d, int col, int row, 00198 char *value); 00199 extern void N_put_array_2d_c_value(N_array_2d * array2d, int col, int row, 00200 CELL value); 00201 extern void N_put_array_2d_f_value(N_array_2d * array2d, int col, int row, 00202 FCELL value); 00203 extern void N_put_array_2d_d_value(N_array_2d * array2d, int col, int row, 00204 DCELL value); 00205 extern int N_is_array_2d_value_null(N_array_2d * array2d, int col, int row); 00206 extern void N_put_array_2d_value_null(N_array_2d * array2d, int col, int row); 00207 extern void N_print_array_2d(N_array_2d * data); 00208 extern void N_print_array_2d_info(N_array_2d * data); 00209 extern void N_copy_array_2d(N_array_2d * source, N_array_2d * target); 00210 extern double N_norm_array_2d(N_array_2d * array1, N_array_2d * array2, 00211 int type); 00212 extern N_array_2d *N_math_array_2d(N_array_2d * array1, N_array_2d * array2, 00213 N_array_2d * result, int type); 00214 extern int N_convert_array_2d_null_to_zero(N_array_2d * a); 00215 extern N_array_2d *N_read_rast_to_array_2d(char *name, N_array_2d * array); 00216 extern void N_write_array_2d_to_rast(N_array_2d * array, char *name); 00217 extern void N_calc_array_2d_stats(N_array_2d * a, double *min, double *max, 00218 double *sum, int *nonzero, int withoffset); 00219 00220 typedef struct 00221 { 00222 int type; /* which raster type FCELL_TYPE, DCELL_TYPE */ 00223 int rows, cols, depths; 00224 int rows_intern, cols_intern, depths_intern; 00225 int offset; /*number of cols/rows/depths offset at each boundary */ 00226 float *fcell_array; /*The data is stored in an one dimensional array internally */ 00227 double *dcell_array; /*The data is stored in an one dimensional array internally */ 00228 } N_array_3d; 00229 00230 extern N_array_3d *N_alloc_array_3d(int cols, int rows, int depths, 00231 int offset, int type); 00232 extern void N_free_array_3d(N_array_3d * data_array); 00233 extern int N_get_array_3d_type(N_array_3d * array3d); 00234 extern void N_get_array_3d_value(N_array_3d * array3d, int col, int row, 00235 int depth, void *value); 00236 extern float N_get_array_3d_f_value(N_array_3d * array3d, int col, int row, 00237 int depth); 00238 extern double N_get_array_3d_d_value(N_array_3d * array3d, int col, int row, 00239 int depth); 00240 extern void N_put_array_3d_value(N_array_3d * array3d, int col, int row, 00241 int depth, char *value); 00242 extern void N_put_array_3d_f_value(N_array_3d * array3d, int col, int row, 00243 int depth, float value); 00244 extern void N_put_array_3d_d_value(N_array_3d * array3d, int col, int row, 00245 int depth, double value); 00246 extern int N_is_array_3d_value_null(N_array_3d * array3d, int col, int row, 00247 int depth); 00248 extern void N_put_array_3d_value_null(N_array_3d * array3d, int col, int row, 00249 int depth); 00250 extern void N_print_array_3d(N_array_3d * data); 00251 extern void N_print_array_3d_info(N_array_3d * data); 00252 extern void N_copy_array_3d(N_array_3d * source, N_array_3d * target); 00253 extern double N_norm_array_3d(N_array_3d * array1, N_array_3d * array2, 00254 int type); 00255 extern N_array_3d *N_math_array_3d(N_array_3d * array1, N_array_3d * array2, 00256 N_array_3d * result, int type); 00257 extern int N_convert_array_3d_null_to_zero(N_array_3d * a); 00258 extern N_array_3d *N_read_rast3d_to_array_3d(char *name, N_array_3d * array, 00259 int mask); 00260 extern void N_write_array_3d_to_rast3d(N_array_3d * array, char *name, 00261 int mask); 00262 extern void N_calc_array_3d_stats(N_array_3d * a, double *min, double *max, 00263 double *sum, int *nonzero, int withoffset); 00264 00265 /* *************************************************************** */ 00266 /* *************** MATRIX ASSEMBLING METHODS ********************* */ 00267 /* *************************************************************** */ 00341 typedef struct 00342 { 00343 int type; 00344 int count; 00345 double C, W, E, N, S, NE, NW, SE, SW, V; 00346 /*top part */ 00347 double T, W_T, E_T, N_T, S_T, NE_T, NW_T, SE_T, SW_T; 00348 /*bottom part */ 00349 double B, W_B, E_B, N_B, S_B, NE_B, NW_B, SE_B, SW_B; 00350 } N_data_star; 00351 00355 typedef struct 00356 { 00357 N_data_star *(*callback) (); 00358 } N_les_callback_3d; 00359 00363 typedef struct 00364 { 00365 N_data_star *(*callback) (); 00366 } N_les_callback_2d; 00367 00368 00369 extern void N_set_les_callback_3d_func(N_les_callback_3d * data, 00370 N_data_star * (*callback_func_3d) ()); 00371 extern void N_set_les_callback_2d_func(N_les_callback_2d * data, 00372 N_data_star * (*callback_func_2d) ()); 00373 extern N_les_callback_3d *N_alloc_les_callback_3d(void); 00374 extern N_les_callback_2d *N_alloc_les_callback_2d(void); 00375 extern N_data_star *N_alloc_5star(void); 00376 extern N_data_star *N_alloc_7star(void); 00377 extern N_data_star *N_alloc_9star(void); 00378 extern N_data_star *N_alloc_27star(void); 00379 extern N_data_star *N_create_5star(double C, double W, double E, double N, 00380 double S, double V); 00381 extern N_data_star *N_create_7star(double C, double W, double E, double N, 00382 double S, double T, double B, double V); 00383 extern N_data_star *N_create_9star(double C, double W, double E, double N, 00384 double S, double NW, double SW, double NE, 00385 double SE, double V); 00386 extern N_data_star *N_create_27star(double C, double W, double E, double N, 00387 double S, double NW, double SW, double NE, 00388 double SE, double T, double W_T, 00389 double E_T, double N_T, double S_T, 00390 double NW_T, double SW_T, double NE_T, 00391 double SE_T, double B, double W_B, 00392 double E_B, double N_B, double S_B, 00393 double NW_B, double SW_B, double NE_B, 00394 double SE_B, double V); 00395 00396 extern N_data_star *N_callback_template_3d(void *data, N_geom_data * geom, 00397 int col, int row, int depth); 00398 extern N_data_star *N_callback_template_2d(void *data, N_geom_data * geom, 00399 int col, int row); 00400 extern N_les *N_assemble_les_3d(int les_type, N_geom_data * geom, 00401 N_array_3d * status, N_array_3d * start_val, 00402 void *data, N_les_callback_3d * callback); 00403 extern N_les *N_assemble_les_3d_active(int les_type, N_geom_data * geom, 00404 N_array_3d * status, 00405 N_array_3d * start_val, void *data, 00406 N_les_callback_3d * callback); 00407 extern N_les *N_assemble_les_3d_dirichlet(int les_type, N_geom_data * geom, 00408 N_array_3d * status, 00409 N_array_3d * start_val, void *data, 00410 N_les_callback_3d * callback); 00411 extern N_les *N_assemble_les_3d_param(int les_type, N_geom_data * geom, 00412 N_array_3d * status, 00413 N_array_3d * start_val, void *data, 00414 N_les_callback_3d * callback, 00415 int cell_type); 00416 extern N_les *N_assemble_les_2d(int les_type, N_geom_data * geom, 00417 N_array_2d * status, N_array_2d * start_val, 00418 void *data, N_les_callback_2d * callback); 00419 extern N_les *N_assemble_les_2d_active(int les_type, N_geom_data * geom, 00420 N_array_2d * status, 00421 N_array_2d * start_val, void *data, 00422 N_les_callback_2d * callback); 00423 extern N_les *N_assemble_les_2d_dirichlet(int les_type, N_geom_data * geom, 00424 N_array_2d * status, 00425 N_array_2d * start_val, void *data, 00426 N_les_callback_2d * callback); 00427 extern N_les *N_assemble_les_2d_param(int les_type, N_geom_data * geom, 00428 N_array_2d * status, 00429 N_array_2d * start_val, void *data, 00430 N_les_callback_2d * callback, 00431 int cell_Type); 00432 00433 extern int N_les_pivot_create(N_les * les); 00434 int N_les_integrate_dirichlet_2d(N_les * les, N_geom_data * geom, 00435 N_array_2d * status, N_array_2d * start_val); 00436 int N_les_integrate_dirichlet_3d(N_les * les, N_geom_data * geom, 00437 N_array_3d * status, N_array_3d * start_val); 00438 00439 /* *************************************************************** */ 00440 /* *************** GPDE STANDARD OPTIONS ************************* */ 00441 /* *************************************************************** */ 00442 00445 typedef enum 00446 { 00447 N_OPT_SOLVER_SYMM, 00448 N_OPT_SOLVER_UNSYMM, 00449 N_OPT_MAX_ITERATIONS, 00450 N_OPT_ITERATION_ERROR, 00451 N_OPT_SOR_VALUE, 00452 N_OPT_CALC_TIME 00453 } N_STD_OPT; 00454 00455 extern struct Option *N_define_standard_option(int opt); 00456 00457 /* *************************************************************** */ 00458 /* *************** GPDE MATHEMATICAL TOOLS *********************** */ 00459 /* *************************************************************** */ 00460 00461 extern double N_calc_arith_mean(double a, double b); 00462 extern double N_calc_arith_mean_n(double *a, int size); 00463 extern double N_calc_geom_mean(double a, double b); 00464 extern double N_calc_geom_mean_n(double *a, int size); 00465 extern double N_calc_harmonic_mean(double a, double b); 00466 extern double N_calc_harmonic_mean_n(double *a, int size); 00467 extern double N_calc_quad_mean(double a, double b); 00468 extern double N_calc_quad_mean_n(double *a, int size); 00469 00470 /* *************************************************************** */ 00471 /* *************** UPWIND STABILIZATION ALGORITHMS *************** */ 00472 /* *************************************************************** */ 00473 00474 extern double N_full_upwinding(double sprod, double distance, double D); 00475 extern double N_exp_upwinding(double sprod, double distance, double D); 00476 00477 00478 /* *************************************************************** */ 00479 /* *************** METHODS FOR GRADIENT CALCULATION ************** */ 00480 /* *************************************************************** */ 00509 typedef struct 00510 { 00511 00512 double NC, SC, WC, EC; 00513 00514 } N_gradient_2d; 00515 00517 typedef struct 00518 { 00519 00520 double NC, SC, WC, EC, TC, BC; 00521 00522 } N_gradient_3d; 00523 00524 00569 typedef struct 00570 { 00571 00572 double NWN, NEN, WC, EC, SWS, SES; 00573 00574 } N_gradient_neighbours_x; 00575 00577 typedef struct 00578 { 00579 00580 double NWW, NEE, NC, SC, SWW, SEE; 00581 00582 } N_gradient_neighbours_y; 00583 00585 typedef struct 00586 { 00587 00588 double NWZ, NZ, NEZ, WZ, CZ, EZ, SWZ, SZ, SEZ; 00589 00590 } N_gradient_neighbours_z; 00591 00593 typedef struct 00594 { 00595 00596 N_gradient_neighbours_x *x; 00597 N_gradient_neighbours_y *y; 00598 00599 } N_gradient_neighbours_2d; 00600 00601 00603 typedef struct 00604 { 00605 00606 N_gradient_neighbours_x *xt; /*top values */ 00607 N_gradient_neighbours_x *xc; /*center values */ 00608 N_gradient_neighbours_x *xb; /*bottom values */ 00609 00610 N_gradient_neighbours_y *yt; /*top values */ 00611 N_gradient_neighbours_y *yc; /*center values */ 00612 N_gradient_neighbours_y *yb; /*bottom values */ 00613 00614 N_gradient_neighbours_z *zt; /*top-center values */ 00615 N_gradient_neighbours_z *zb; /*bottom-center values */ 00616 00617 } N_gradient_neighbours_3d; 00618 00619 00621 typedef struct 00622 { 00623 00624 N_array_2d *x_array; 00625 N_array_2d *y_array; 00626 int cols, rows; 00627 double min, max, mean, sum; 00628 int nonull; 00629 00630 } N_gradient_field_2d; 00631 00633 typedef struct 00634 { 00635 00636 N_array_3d *x_array; 00637 N_array_3d *y_array; 00638 N_array_3d *z_array; 00639 int cols, rows, depths; 00640 double min, max, mean, sum; 00641 int nonull; 00642 00643 } N_gradient_field_3d; 00644 00645 00646 extern N_gradient_2d *N_alloc_gradient_2d(void); 00647 extern void N_free_gradient_2d(N_gradient_2d * grad); 00648 extern N_gradient_2d *N_create_gradient_2d(double NC, double SC, double WC, 00649 double EC); 00650 extern int N_copy_gradient_2d(N_gradient_2d * source, N_gradient_2d * target); 00651 extern N_gradient_2d *N_get_gradient_2d(N_gradient_field_2d * field, 00652 N_gradient_2d * gradient, int col, 00653 int row); 00654 00655 extern N_gradient_3d *N_alloc_gradient_3d(void); 00656 extern void N_free_gradient_3d(N_gradient_3d * grad); 00657 extern N_gradient_3d *N_create_gradient_3d(double NC, double SC, double WC, 00658 double EC, double TC, double BC); 00659 extern int N_copy_gradient_3d(N_gradient_3d * source, N_gradient_3d * target); 00660 extern N_gradient_3d *N_get_gradient_3d(N_gradient_field_3d * field, 00661 N_gradient_3d * gradient, int col, 00662 int row, int depth); 00663 00664 extern N_gradient_neighbours_x *N_alloc_gradient_neighbours_x(void); 00665 extern void N_free_gradient_neighbours_x(N_gradient_neighbours_x * grad); 00666 extern N_gradient_neighbours_x *N_create_gradient_neighbours_x(double NWN, 00667 double NEN, 00668 double WC, 00669 double EC, 00670 double SWS, 00671 double SES); 00672 extern int N_copy_gradient_neighbours_x(N_gradient_neighbours_x * source, 00673 N_gradient_neighbours_x * target); 00674 00675 extern N_gradient_neighbours_y *N_alloc_gradient_neighbours_y(void); 00676 extern void N_free_gradient_neighbours_y(N_gradient_neighbours_y * grad); 00677 extern N_gradient_neighbours_y *N_create_gradient_neighbours_y(double NWW, 00678 double NEE, 00679 double NC, 00680 double SC, 00681 double SWW, 00682 double SEE); 00683 extern int N_copy_gradient_neighbours_y(N_gradient_neighbours_y * source, 00684 N_gradient_neighbours_y * target); 00685 00686 extern N_gradient_neighbours_z *N_alloc_gradient_neighbours_z(void); 00687 extern void N_free_gradient_neighbours_z(N_gradient_neighbours_z * grad); 00688 extern N_gradient_neighbours_z *N_create_gradient_neighbours_z(double NWZ, 00689 double NZ, 00690 double NEZ, 00691 double WZ, 00692 double CZ, 00693 double EZ, 00694 double SWZ, 00695 double SZ, 00696 double SEZ); 00697 extern int N_copy_gradient_neighbours_z(N_gradient_neighbours_z * source, 00698 N_gradient_neighbours_z * target); 00699 00700 extern N_gradient_neighbours_2d *N_alloc_gradient_neighbours_2d(void); 00701 extern void N_free_gradient_neighbours_2d(N_gradient_neighbours_2d * grad); 00702 extern N_gradient_neighbours_2d 00703 * N_create_gradient_neighbours_2d(N_gradient_neighbours_x * x, 00704 N_gradient_neighbours_y * y); 00705 extern int N_copy_gradient_neighbours_2d(N_gradient_neighbours_2d * source, 00706 N_gradient_neighbours_2d * target); 00707 extern N_gradient_neighbours_2d 00708 * N_get_gradient_neighbours_2d(N_gradient_field_2d * field, 00709 N_gradient_neighbours_2d * gradient, 00710 int col, int row); 00711 00712 00713 extern N_gradient_neighbours_3d *N_alloc_gradient_neighbours_3d(void); 00714 extern void N_free_gradient_neighbours_3d(N_gradient_neighbours_3d * grad); 00715 extern N_gradient_neighbours_3d 00716 * N_create_gradient_neighbours_3d(N_gradient_neighbours_x * xt, 00717 N_gradient_neighbours_x * xc, 00718 N_gradient_neighbours_x * xb, 00719 N_gradient_neighbours_y * yt, 00720 N_gradient_neighbours_y * yc, 00721 N_gradient_neighbours_y * yb, 00722 N_gradient_neighbours_z * zt, 00723 N_gradient_neighbours_z * zb); 00724 extern int N_copy_gradient_neighbours_3d(N_gradient_neighbours_3d * source, 00725 N_gradient_neighbours_3d * target); 00726 00727 extern void N_print_gradient_field_2d_info(N_gradient_field_2d * field); 00728 extern void N_calc_gradient_field_2d_stats(N_gradient_field_2d * field); 00729 00730 00731 extern N_gradient_field_2d *N_alloc_gradient_field_2d(int cols, int rows); 00732 extern void N_free_gradient_field_2d(N_gradient_field_2d * field); 00733 extern int N_copy_gradient_field_2d(N_gradient_field_2d * source, 00734 N_gradient_field_2d * target); 00735 extern N_gradient_field_2d *N_compute_gradient_field_2d(N_array_2d * pot, 00736 N_array_2d * weight_x, 00737 N_array_2d * weight_y, 00738 N_geom_data * geom, 00739 N_gradient_field_2d * 00740 gradfield); 00741 extern void N_compute_gradient_field_components_2d(N_gradient_field_2d * 00742 field, N_array_2d * x_comp, 00743 N_array_2d * y_comp); 00744 00745 extern void N_print_gradient_field_3d_info(N_gradient_field_3d * field); 00746 extern void N_calc_gradient_field_3d_stats(N_gradient_field_3d * field); 00747 00748 extern N_gradient_field_3d *N_alloc_gradient_field_3d(int cols, int rows, 00749 int depths); 00750 extern void N_free_gradient_field_3d(N_gradient_field_3d * field); 00751 extern int N_copy_gradient_field_3d(N_gradient_field_3d * source, 00752 N_gradient_field_3d * target); 00753 extern N_gradient_field_3d *N_compute_gradient_field_3d(N_array_3d * pot, 00754 N_array_3d * weight_x, 00755 N_array_3d * weight_y, 00756 N_array_3d * weight_z, 00757 N_geom_data * geom, 00758 N_gradient_field_3d * 00759 gradfield); 00760 extern void N_compute_gradient_field_components_3d(N_gradient_field_3d * 00761 field, N_array_3d * x_comp, 00762 N_array_3d * y_comp, 00763 N_array_3d * z_comp); 00764 00765 #endif