rm files
parent
6ac4de81a9
commit
562b948508
@ -1,68 +0,0 @@
|
|||||||
#include <stdio.h>
|
|
||||||
#include <stddef.h>
|
|
||||||
#include <string.h>
|
|
||||||
#include <stdarg.h>
|
|
||||||
#include <stdlib.h>
|
|
||||||
#include <math.h>
|
|
||||||
#include "geostat.h"
|
|
||||||
#include "pressure.h"
|
|
||||||
|
|
||||||
/* Build_pressure */
|
|
||||||
/* build pressure in spectral domain */
|
|
||||||
|
|
||||||
void build_pressure(int n[3],struct grid_mod grid,struct pressure_mod gradient,double *realization,double *ireal,double *pressure,double *ipressure)
|
|
||||||
{
|
|
||||||
int i,j,k,maille1;
|
|
||||||
double ki,kj,kk;
|
|
||||||
double coeff,temp,temp2;
|
|
||||||
|
|
||||||
printf("build pdsdsressure pdslis\n");
|
|
||||||
|
|
||||||
/* pressure calculation in the spectral domain*/
|
|
||||||
for ( k = 1; k <= n[2]; k++) {
|
|
||||||
for (j = 1; j <= n[1]; j++) {
|
|
||||||
for (i = 1; i <= n[0]; i++) {
|
|
||||||
maille1 = i+(j-1+(k-1)*n[1])*n[0];
|
|
||||||
if (n[0]==1)
|
|
||||||
{
|
|
||||||
ki=0.;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
ki =(double)(i-1)/(double)(grid.NX*grid.DX);
|
|
||||||
}
|
|
||||||
if (n[1]==1)
|
|
||||||
{
|
|
||||||
kj=0.;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
kj =(double)j/(double)(grid.NY*grid.DY);
|
|
||||||
}
|
|
||||||
if (n[2]==1)
|
|
||||||
{
|
|
||||||
kk=0.;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
kk =(double)k/(double)(grid.NZ*grid.DZ);
|
|
||||||
}
|
|
||||||
|
|
||||||
coeff = (gradient.x*ki+gradient.y*kj+gradient.z*kk)/(ki*ki+kj*kj+kk*kk);
|
|
||||||
temp = realization[maille1];
|
|
||||||
temp2= ireal[maille1];
|
|
||||||
if (maille1==1)
|
|
||||||
{
|
|
||||||
pressure[maille1] =0.;
|
|
||||||
ipressure[maille1]=0.;
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
pressure[maille1] =-1./(2*3.14)*coeff*temp2;
|
|
||||||
ipressure[maille1]= 1./(2*3.14)*coeff*temp;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
return;
|
|
||||||
}
|
|
@ -1,119 +0,0 @@
|
|||||||
#include <stdio.h>
|
|
||||||
#include <stddef.h>
|
|
||||||
#include <stdarg.h>
|
|
||||||
#include <stdlib.h>
|
|
||||||
#include <string.h>
|
|
||||||
#include <math.h>
|
|
||||||
#include "geostat.h"
|
|
||||||
#include "pressure.h"
|
|
||||||
|
|
||||||
/* Build_velocity */
|
|
||||||
/* Build velocity in spectral domain */
|
|
||||||
|
|
||||||
void build_velocity(int n[3],struct grid_mod grid,struct statistic_mod stat,struct pressure_mod gradient,double *realization,double *ireal,double *xvelocity,double *ixvelocity,int direction)
|
|
||||||
{
|
|
||||||
int i,j,k,maille1;
|
|
||||||
int x,y,z;
|
|
||||||
double ki,kj,kk;
|
|
||||||
double coeff,temp,temp2;
|
|
||||||
double grad;
|
|
||||||
|
|
||||||
printf("build velocity\n");
|
|
||||||
|
|
||||||
/* x-direction velocity calculation in the spectral domain*/
|
|
||||||
switch(direction)
|
|
||||||
{
|
|
||||||
case 1:
|
|
||||||
grad = gradient.x;
|
|
||||||
x=1;
|
|
||||||
y=0;
|
|
||||||
z=0;
|
|
||||||
break;
|
|
||||||
case 2:
|
|
||||||
grad = gradient.y;
|
|
||||||
x=0;
|
|
||||||
y=1;
|
|
||||||
z=0;
|
|
||||||
break;
|
|
||||||
case 3:
|
|
||||||
grad = gradient.z;
|
|
||||||
x=0;
|
|
||||||
y=0;
|
|
||||||
z=1;
|
|
||||||
break;
|
|
||||||
default :
|
|
||||||
printf("build_velocity.c: wrong velocity direction!!! direction: %d\n",direction);
|
|
||||||
break;
|
|
||||||
}
|
|
||||||
|
|
||||||
for ( k = 1; k <= n[2]; k++) {
|
|
||||||
for (j = 1; j <= n[1]; j++) {
|
|
||||||
for (i = 1; i <= n[0]; i++) {
|
|
||||||
maille1 = i+(j-1+(k-1)*n[1])*n[0];
|
|
||||||
/* if (i<n[0]/2) { */
|
|
||||||
/* ki =(double)i/(double)(grid.NX*grid.DX); */
|
|
||||||
/* } else { */
|
|
||||||
/* ki =-(double)(n[0]-i+1)/(double)(grid.NX*grid.DX); */
|
|
||||||
/* } */
|
|
||||||
/* if (j<n[1]/2) { */
|
|
||||||
/* kj =(double)j/(double)(grid.NY*grid.DY); */
|
|
||||||
/* } else { */
|
|
||||||
/* kj =-(double)(n[1]-j+1)/(double)(grid.NY*grid.DY); */
|
|
||||||
/* } */
|
|
||||||
/* if (k<n[2]/2) { */
|
|
||||||
/* kk =(double)k/(double)(grid.NZ*grid.DZ); */
|
|
||||||
/* } else { */
|
|
||||||
/* kk =-(double)(n[2]-k+1)/(double)(grid.NZ*grid.DZ); */
|
|
||||||
/* } */
|
|
||||||
|
|
||||||
if (n[0]==1)
|
|
||||||
{
|
|
||||||
ki=0;
|
|
||||||
}
|
|
||||||
else if (i<n[0]/2)
|
|
||||||
{
|
|
||||||
ki =(double)i/(double)(grid.NX*grid.DX);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
ki =-(double)(n[0]-i+1)/(double)(grid.NX*grid.DX);
|
|
||||||
}
|
|
||||||
|
|
||||||
if (n[1]==1)
|
|
||||||
{
|
|
||||||
kj=0;
|
|
||||||
}
|
|
||||||
else if (j<n[1]/2)
|
|
||||||
{
|
|
||||||
kj =(double)j/(double)(grid.NY*grid.DY);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
kj =-(double)(n[1]-j+1)/(double)(grid.NY*grid.DY);
|
|
||||||
}
|
|
||||||
|
|
||||||
if (n[2]==1)
|
|
||||||
{
|
|
||||||
kk=0;
|
|
||||||
}
|
|
||||||
else if (k<n[2]/2)
|
|
||||||
{
|
|
||||||
kk =(double)k/(double)(grid.NZ*grid.DZ);
|
|
||||||
}
|
|
||||||
else
|
|
||||||
{
|
|
||||||
kk =-(double)(n[2]-k+1)/(double)(grid.NZ*grid.DZ);
|
|
||||||
}
|
|
||||||
|
|
||||||
coeff = (gradient.x*ki+gradient.y*kj+gradient.z*kk)/(ki*ki+kj*kj+kk*kk);
|
|
||||||
temp = realization[maille1];
|
|
||||||
temp2= ireal[maille1];
|
|
||||||
|
|
||||||
xvelocity[maille1] = stat.mean[0]*(grad-(ki*x+kj*y+kk*z)*coeff)*temp;
|
|
||||||
ixvelocity[maille1]= stat.mean[0]*(grad-(ki*x+kj*y+kk*z)*coeff)*temp2;
|
|
||||||
}
|
|
||||||
}
|
|
||||||
}
|
|
||||||
|
|
||||||
return;
|
|
||||||
}
|
|
@ -1,33 +0,0 @@
|
|||||||
/*calculates C.x/
|
|
||||||
/* C : symetric positive-definite matrix recorded */
|
|
||||||
/* (per raws) as a vector with only components */
|
|
||||||
/* Cij so that j <= i, i = [0...n-1]*/
|
|
||||||
/* x : vector, xi avec i = [0...n-1]*/
|
|
||||||
/* b : vector, bi avec i = [0...n-1]*/
|
|
||||||
/* n : dimension of matrix Cij*/
|
|
||||||
/* */
|
|
||||||
/* The solution vector is returned in b*/
|
|
||||||
|
|
||||||
void mat_vec(double *C, double *x, double *b, int n)
|
|
||||||
|
|
||||||
{
|
|
||||||
int i,k,l;
|
|
||||||
double temp;
|
|
||||||
|
|
||||||
for (i = 0; i < n; i++) {
|
|
||||||
temp = 0.;
|
|
||||||
|
|
||||||
for (k = 0; k < n; k++) {
|
|
||||||
if ( k <= i ) {
|
|
||||||
l = k+i*(i+1)/2;
|
|
||||||
} else {
|
|
||||||
l = i+k*(k+1)/2;
|
|
||||||
}
|
|
||||||
temp += C[l]*x[k];
|
|
||||||
}
|
|
||||||
|
|
||||||
b[i] = temp;
|
|
||||||
}
|
|
||||||
|
|
||||||
return;
|
|
||||||
}
|
|
Loading…
Reference in New Issue