You cannot select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
100 lines
3.5 KiB
C
100 lines
3.5 KiB
C
#include "geostat.h"
|
|
#include "log.h"
|
|
#include <time.h>
|
|
|
|
/*builds the sampled covariance function*/
|
|
/*dimensions are even*/
|
|
void covariance(double* covar, struct vario_mod variogram, struct grid_mod mesh, int n[3]) {
|
|
clock_t t = clock();
|
|
|
|
int i, j, k, maille, n2[3], symmetric;
|
|
double di, dj, dk;
|
|
|
|
log_info("RESULT = in progress");
|
|
|
|
for (i = 0; i < 3; i++)
|
|
n2[i] = n[i] / 2;
|
|
|
|
for (i = 0; i <= n2[0]; i++) {
|
|
for (j = 0; j <= n2[1]; j++) {
|
|
for (k = 0; k <= n2[2]; k++) {
|
|
|
|
/*area 1*/
|
|
maille = 1 + i + n[0] * (j + n[1] * k);
|
|
di = (double)i * mesh.DX;
|
|
dj = (double)j * mesh.DY;
|
|
dk = (double)k * mesh.DZ;
|
|
covar[maille] = (double)cov_value(variogram, di, dj, dk);
|
|
|
|
if (k > 0 && k < n2[2] && j > 0 && j < n2[1] && i > 0 && i < n2[0]) {
|
|
/*area 2*/
|
|
symmetric = 1 + n[0] - i + n[0] * (n[1] - j + n[1] * (n[2] - k));
|
|
covar[symmetric] = covar[maille];
|
|
}
|
|
|
|
if (i > 0 && i < n2[0]) {
|
|
/*area 4*/
|
|
di = -(double)i * mesh.DX;
|
|
dj = (double)j * mesh.DY;
|
|
dk = (double)k * mesh.DZ;
|
|
maille = 1 + (n[0] - i) + n[0] * (j + n[1] * k);
|
|
covar[maille] = (double)cov_value(variogram, di, dj, dk);
|
|
}
|
|
|
|
if (k > 0 && k < n2[2] && j > 0 && j < n2[1]) {
|
|
/*area 8*/
|
|
symmetric = 1 + i + n[0] * (n[1] - j + n[1] * (n[2] - k));
|
|
covar[symmetric] = covar[maille];
|
|
}
|
|
|
|
if (i > 0 && i < n2[0] && j > 0 && j < n2[1]) {
|
|
/*area 5*/
|
|
di = -(double)i * mesh.DX;
|
|
dj = -(double)j * mesh.DY;
|
|
dk = (double)k * mesh.DZ;
|
|
maille = 1 + (n[0] - i) + n[0] * (n[1] - j + n[1] * k);
|
|
covar[maille] = (double)cov_value(variogram, di, dj, dk);
|
|
}
|
|
|
|
if (k > 0 && k < n2[2]) {
|
|
/*area 6*/
|
|
symmetric = 1 + i + n[0] * (j + n[1] * (n[2] - k));
|
|
covar[symmetric] = covar[maille];
|
|
}
|
|
|
|
if (j > 0 && j < n2[1]) {
|
|
/*area 3*/
|
|
di = (double)i * mesh.DX;
|
|
dj = -(double)j * mesh.DY;
|
|
dk = (double)k * mesh.DZ;
|
|
maille = 1 + i + n[0] * (n[1] - j + n[1] * k);
|
|
covar[maille] = (double)cov_value(variogram, di, dj, dk);
|
|
}
|
|
|
|
if (k > 0 && k < n2[2] && i > 0 && i < n2[0]) {
|
|
/*area 7*/
|
|
symmetric = 1 + n[0] - i + n[0] * (j + n[1] * (n[2] - k));
|
|
covar[symmetric] = covar[maille];
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
t = clock() - t;
|
|
double time_taken = ((double)t)/CLOCKS_PER_SEC; // calculate the elapsed time
|
|
|
|
double* total_ram = malloc(sizeof(double));
|
|
getTotalVirtualMem(total_ram);
|
|
|
|
double* used_ram = malloc(sizeof(double));
|
|
getVirtualMemUsed(used_ram);
|
|
|
|
log_info("TOTAL VIRTUAL MEM = %5.1f MB, USED VIRTUAL MEM = %5.1f MB, USED VIRTUAL MEM BY CURRENT PROCESS = %d MB",
|
|
*total_ram, *used_ram, getVirtualMemUsedByCurrentProcess());
|
|
|
|
free(total_ram);
|
|
free(used_ram);
|
|
|
|
log_info("RESULT = success, di = %f, dj = %f, dk = %f, ELAPSED = %f seconds", di, dj, dk, time_taken);
|
|
}
|