GRASS Programmer's Manual
6.4.2(2012)
Main Page
Related Pages
Namespaces
Data Structures
Files
File List
Globals
All
Data Structures
Namespaces
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
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
);
214
extern
int
N_convert_array_2d_null_to_zero
(
N_array_2d
* a);
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
);
257
extern
int
N_convert_array_3d_null_to_zero
(
N_array_3d
* a);
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) ();
358
}
N_les_callback_3d
;
359
363
typedef
struct
364
{
365
N_data_star
*(*callback) ();
366
}
N_les_callback_2d
;
367
368
369
extern
void
N_set_les_callback_3d_func
(
N_les_callback_3d
*
data
,
370
N_data_star
* (*callback_func_3d) ());
371
extern
void
N_set_les_callback_2d_func
(
N_les_callback_2d
* data,
372
N_data_star
* (*callback_func_2d) ());
373
extern
N_les_callback_3d
*
N_alloc_les_callback_3d
(
void
);
374
extern
N_les_callback_2d
*
N_alloc_les_callback_2d
(
void
);
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);
434
int
N_les_integrate_dirichlet_2d
(
N_les
* les,
N_geom_data
* geom,
435
N_array_2d
* status,
N_array_2d
* start_val);
436
int
N_les_integrate_dirichlet_3d
(
N_les
* les,
N_geom_data
* geom,
437
N_array_3d
* status,
N_array_3d
* start_val);
438
439
/* *************************************************************** */
440
/* *************** GPDE STANDARD OPTIONS ************************* */
441
/* *************************************************************** */
442
445
typedef
enum
446
{
447
N_OPT_SOLVER_SYMM
,
448
N_OPT_SOLVER_UNSYMM
,
449
N_OPT_MAX_ITERATIONS
,
450
N_OPT_ITERATION_ERROR
,
451
N_OPT_SOR_VALUE
,
452
N_OPT_CALC_TIME
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
574
}
N_gradient_neighbours_x
;
575
577
typedef
struct
578
{
579
580
double
NWW, NEE, NC,
SC
,
SWW
, SEE;
581
582
}
N_gradient_neighbours_y
;
583
585
typedef
struct
586
{
587
588
double
NWZ, NZ, NEZ,
WZ
, CZ, EZ, SWZ, SZ, SEZ;
589
590
}
N_gradient_neighbours_z
;
591
593
typedef
struct
594
{
595
596
N_gradient_neighbours_x
*
x
;
597
N_gradient_neighbours_y
*
y
;
598
599
}
N_gradient_neighbours_2d
;
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
617
}
N_gradient_neighbours_3d
;
618
619
621
typedef
struct
622
{
623
624
N_array_2d
*
x_array
;
625
N_array_2d
*
y_array
;
626
int
cols
,
rows
;
627
double
min
,
max
, mean,
sum
;
628
int
nonull
;
629
630
}
N_gradient_field_2d
;
631
633
typedef
struct
634
{
635
636
N_array_3d
*
x_array
;
637
N_array_3d
*
y_array
;
638
N_array_3d
*
z_array
;
639
int
cols
,
rows
, depths;
640
double
min
,
max
, mean,
sum
;
641
int
nonull
;
642
643
}
N_gradient_field_3d
;
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);
651
extern
N_gradient_2d
*
N_get_gradient_2d
(
N_gradient_field_2d
* field,
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);
660
extern
N_gradient_3d
*
N_get_gradient_3d
(
N_gradient_field_3d
* field,
661
N_gradient_3d
* gradient,
int
col,
662
int
row,
int
depth);
663
664
extern
N_gradient_neighbours_x
*
N_alloc_gradient_neighbours_x
(
void
);
665
extern
void
N_free_gradient_neighbours_x
(
N_gradient_neighbours_x
* grad);
666
extern
N_gradient_neighbours_x
*
N_create_gradient_neighbours_x
(
double
NWN,
667
double
NEN,
668
double
WC,
669
double
EC,
670
double
SWS,
671
double
SES);
672
extern
int
N_copy_gradient_neighbours_x
(
N_gradient_neighbours_x
* source,
673
N_gradient_neighbours_x
* target);
674
675
extern
N_gradient_neighbours_y
*
N_alloc_gradient_neighbours_y
(
void
);
676
extern
void
N_free_gradient_neighbours_y
(
N_gradient_neighbours_y
* grad);
677
extern
N_gradient_neighbours_y
*
N_create_gradient_neighbours_y
(
double
NWW,
678
double
NEE,
679
double
NC,
680
double
SC
,
681
double
SWW,
682
double
SEE);
683
extern
int
N_copy_gradient_neighbours_y
(
N_gradient_neighbours_y
* source,
684
N_gradient_neighbours_y
* target);
685
686
extern
N_gradient_neighbours_z
*
N_alloc_gradient_neighbours_z
(
void
);
687
extern
void
N_free_gradient_neighbours_z
(
N_gradient_neighbours_z
* grad);
688
extern
N_gradient_neighbours_z
*
N_create_gradient_neighbours_z
(
double
NWZ,
689
double
NZ,
690
double
NEZ,
691
double
WZ,
692
double
CZ,
693
double
EZ,
694
double
SWZ,
695
double
SZ,
696
double
SEZ);
697
extern
int
N_copy_gradient_neighbours_z
(
N_gradient_neighbours_z
* source,
698
N_gradient_neighbours_z
* target);
699
700
extern
N_gradient_neighbours_2d
*
N_alloc_gradient_neighbours_2d
(
void
);
701
extern
void
N_free_gradient_neighbours_2d
(
N_gradient_neighbours_2d
* grad);
702
extern
N_gradient_neighbours_2d
703
*
N_create_gradient_neighbours_2d
(
N_gradient_neighbours_x
* x,
704
N_gradient_neighbours_y
* y);
705
extern
int
N_copy_gradient_neighbours_2d
(
N_gradient_neighbours_2d
* source,
706
N_gradient_neighbours_2d
* target);
707
extern
N_gradient_neighbours_2d
708
*
N_get_gradient_neighbours_2d
(
N_gradient_field_2d
* field,
709
N_gradient_neighbours_2d
* gradient,
710
int
col,
int
row);
711
712
713
extern
N_gradient_neighbours_3d
*
N_alloc_gradient_neighbours_3d
(
void
);
714
extern
void
N_free_gradient_neighbours_3d
(
N_gradient_neighbours_3d
* grad);
715
extern
N_gradient_neighbours_3d
716
*
N_create_gradient_neighbours_3d
(
N_gradient_neighbours_x
* xt,
717
N_gradient_neighbours_x
* xc,
718
N_gradient_neighbours_x
* xb,
719
N_gradient_neighbours_y
* yt,
720
N_gradient_neighbours_y
* yc,
721
N_gradient_neighbours_y
* yb,
722
N_gradient_neighbours_z
* zt,
723
N_gradient_neighbours_z
* zb);
724
extern
int
N_copy_gradient_neighbours_3d
(
N_gradient_neighbours_3d
* source,
725
N_gradient_neighbours_3d
* target);
726
727
extern
void
N_print_gradient_field_2d_info
(
N_gradient_field_2d
* field);
728
extern
void
N_calc_gradient_field_2d_stats
(
N_gradient_field_2d
* field);
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);
733
extern
int
N_copy_gradient_field_2d
(
N_gradient_field_2d
* source,
734
N_gradient_field_2d
* target);
735
extern
N_gradient_field_2d
*
N_compute_gradient_field_2d
(
N_array_2d
* pot,
736
N_array_2d
* weight_x,
737
N_array_2d
* weight_y,
738
N_geom_data
* geom,
739
N_gradient_field_2d
*
740
gradfield);
741
extern
void
N_compute_gradient_field_components_2d
(
N_gradient_field_2d
*
742
field,
N_array_2d
* x_comp,
743
N_array_2d
* y_comp);
744
745
extern
void
N_print_gradient_field_3d_info
(
N_gradient_field_3d
* field);
746
extern
void
N_calc_gradient_field_3d_stats
(
N_gradient_field_3d
* field);
747
748
extern
N_gradient_field_3d
*
N_alloc_gradient_field_3d
(
int
cols
,
int
rows,
749
int
depths);
750
extern
void
N_free_gradient_field_3d
(
N_gradient_field_3d
* field);
751
extern
int
N_copy_gradient_field_3d
(
N_gradient_field_3d
* source,
752
N_gradient_field_3d
* target);
753
extern
N_gradient_field_3d
*
N_compute_gradient_field_3d
(
N_array_3d
* pot,
754
N_array_3d
* weight_x,
755
N_array_3d
* weight_y,
756
N_array_3d
* weight_z,
757
N_geom_data
* geom,
758
N_gradient_field_3d
*
759
gradfield);
760
extern
void
N_compute_gradient_field_components_3d
(
N_gradient_field_3d
*
761
field,
N_array_3d
* x_comp,
762
N_array_3d
* y_comp,
763
N_array_3d
* z_comp);
764
765
#endif
lib
gpde
N_pde.h
Generated on Wed Jun 6 2012 14:04:26 for GRASS Programmer's Manual by
1.8.1