sandbox/yixian/silo.c
# we solve a case of 2D silo with orifice at the bottom by using the poisson equation \displaystyle \nabla^2 c = 0 with a flux at the top of the silo This system can be solved with the poisson solver.
#include "grid/multigrid.h"
#include "poisson.h"
#include "utils.h"
Concentration source, and some obvious variables
scalar c[],s[];
double dx,dy,W;
Boundary conditions: a flux at the top of silo an orifice with size W at the bottom of the silo
c[top] = neumann(1);
c[bottom] = fabs(x)<= W/2. ? dirichlet(0): neumann(0);
c[right] = neumann(0);
c[left] = neumann(0);
## Parameters The size of the domain L0
. 0<x<L0 0<y<L0
int main() {
L0 = 1.;
X0 = -L0/2.;
Y0 = -L0/2.;
N = 64;
dx = L0/N;
dy = L0/N;
W = 0.25;
init_grid (N);
##Solve Laplace Initialisation of the solution of laplace equation with poisson
\displaystyle
\nabla^2 c = s
with a source term s=0
Save the results
numerical solution
FILE * fpc = fopen("Fcsilo.txt", "w");
for (double x =-L0/2 ; x <= L0/2; x += dx){
for (double y = -L0/2 ; y <= L0/2; y += dy){
fprintf (fpc, "%g %g %g \n",
x, y, interpolate (c, x, y));}
fprintf (fpc,"\n");}
fclose(fpc);
fprintf(stdout," end\n");
}
##Run Then compile and run:
bash qcc -O2 -Wall -o silo silo.c -lm; ./silo ##Results
set hidden3d
sp 'Fcsilo.txt'
reset
set pm3d map
set size square
set palette rgbformulae 22,13,-31
sp 'Fcsilo.txt'
##Bibliography
Version 1: 4 december 2015/ 14:58