HyPar  1.0
Finite-Difference Hyperbolic-Parabolic PDE Solver on Cartesian Grids
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
LinearADRAdvectionField.c File Reference

Contains the function to read in a non-constant advection field. More...

#include <stdlib.h>
#include <math.h>
#include <string.h>
#include <basic.h>
#include <common.h>
#include <arrayfunctions.h>
#include <io.h>
#include <physicalmodels/linearadr.h>
#include <hypar.h>
#include <mpivars.h>

Go to the source code of this file.

Functions

int LinearADRAdvectionField (void *s, void *m, int idx, int nsims, int *dim_adv)
 

Detailed Description

Contains the function to read in a non-constant advection field.

Author
Debojyoti Ghosh

Definition in file LinearADRAdvectionField.c.

Function Documentation

int LinearADRAdvectionField ( void *  s,
void *  m,
int  idx,
int  nsims,
int *  dim_adv 
)

Set the advection field over the domain - reads the spatially-varying advection data from a file, if available. The array to store the field must already be allocated.

For a single simulation, it reads in the data from a file. The data must have the same grid dimensions as the solution.

For an ensemble simulation, it reads in the advection field from the file corresponding to the idx (index of this simulation). The data in this file must have the same grid dimensions as the domain it is being read in for.

For a sparse grids simulation, it reads in the advection field from the file with data that has the grid dimensions as the full grid. The field on the current grid is obtained by interpolation.

Parameters
sSolver object of type HyPar
mMPI object of type MPIVariables
idxIndex of this simulation
nsimsTotal number of simulations
dim_advGrid dimensions of the advection field stored in file

Definition at line 33 of file LinearADRAdvectionField.c.

39 {
40  HyPar *solver = (HyPar*) s;
41  MPIVariables *mpi = (MPIVariables*) m;
42  LinearADR *param = (LinearADR*) solver->physics;
43  double *adv = param->a;
44 
45  /* sanity checks */
46  if (param->constant_advection == 1) {
47  fprintf(stderr,"Error in LinearADRAdvectionField(): param->constant_advection == 1!\n");
48  return(1);
49  }
50  if (!strcmp(param->adv_filename,"none")) {
51  fprintf(stderr,"Error in LinearADRAdvectionField(): param->adv_filename is \"none\"!\n");
52  return(1);
53  }
54 
55  int d, done, flag;
56  int *dim = solver->dim_local;
57  int ghosts = solver->ghosts;
58 
59  /* number of advection field components at a grid point is
60  * number of spatial dimensions times number of solution
61  * components */
62  int adv_nvar = solver->ndims * solver->nvars;
63  if (param->adv_arr_size != adv_nvar*solver->npoints_local_wghosts) {
64  fprintf(stderr,"Error in LinearADRAdvectionField(): Incorrect adv_arr_size\n");
65  return(1);
66  }
67 
68  char fname_root[_MAX_STRING_SIZE_];
69  strcpy(fname_root, param->adv_filename);
70  if (idx >= 0) {
71  if (nsims > 1) {
72  char index[_MAX_STRING_SIZE_];
73  GetStringFromInteger(idx, index, (int)log10(nsims)+1);
74  strcat(fname_root, "_");
75  strcat(fname_root, index);
76  }
77  }
78 
79  /* read spatially-varying advection field from provided file */
80  if (dim_adv == NULL) {
81  int ierr = ReadArray( solver->ndims,
82  adv_nvar,
83  solver->dim_global,
84  solver->dim_local,
85  solver->ghosts,
86  solver,
87  mpi,
88  NULL,
89  adv,
90  fname_root,
91  &flag);
92  if (ierr) {
93  fprintf(stderr,"Error in LinearADRAdvectionField()\n");
94  fprintf(stderr," ReadArray() returned error!\n");
95  return ierr;
96  }
97  } else {
98  int ierr = ReadArraywInterp(solver->ndims,
99  adv_nvar,
100  solver->dim_global,
101  solver->dim_local,
102  dim_adv,
103  solver->ghosts,
104  solver,
105  mpi,
106  NULL,
107  adv,
108  fname_root,
109  &flag);
110  if (ierr) {
111  fprintf(stderr,"Error in LinearADRAdvectionField()\n");
112  fprintf(stderr," ReadArraywInterp() returned error!\n");
113  return ierr;
114  }
115  }
116 
117  if (!flag) {
118  /* if advection file not available, set it to zero */
119  _ArraySetValue_(adv,solver->npoints_local_wghosts*adv_nvar,0.0);
120  }
121 
122  /* if parallel, exchange MPI-boundary ghost point data */
124  adv_nvar,
125  solver->dim_local,
126  solver->ghosts,
127  mpi,
128  adv );
129 
130  /* Along each dimension:
131  - If boundaries are periodic, then set ghost points accordingly if
132  number of MPI ranks along that dimension is 1 (if more than 1,
133  MPIExchangeBoundariesnD() called above has taken care of it).
134  - Else, extrapolate at the boundaries. */
135  int indexb[solver->ndims],
136  indexi[solver->ndims],
137  bounds[solver->ndims],
138  offset[solver->ndims];
139  for (d = 0; d < solver->ndims; d++) {
140  if (solver->isPeriodic[d] && (mpi->iproc[d] == 1)) {
141  _ArrayCopy1D_(dim,bounds,solver->ndims); bounds[d] = ghosts;
142  /* left boundary */
143  done = 0; _ArraySetValue_(indexb,solver->ndims,0);
144  _ArraySetValue_(offset,solver->ndims,0); offset[d] = -ghosts;
145  while (!done) {
146  _ArrayCopy1D_(indexb,indexi,solver->ndims); indexi[d] = indexb[d] + dim[d] - ghosts;
147  int p1; _ArrayIndex1DWO_(solver->ndims,dim,indexb,offset,ghosts,p1);
148  int p2; _ArrayIndex1D_ (solver->ndims,dim,indexi,ghosts,p2);
149  _ArrayCopy1D_((adv+p2*adv_nvar), (adv+p1*adv_nvar), adv_nvar);
150  _ArrayIncrementIndex_(solver->ndims,bounds,indexb,done);
151  }
152  /* right boundary */
153  done = 0; _ArraySetValue_(indexb,solver->ndims,0);
154  _ArraySetValue_(offset,solver->ndims,0); offset[d] = dim[d];
155  while (!done) {
156  _ArrayCopy1D_(indexb,indexi,solver->ndims);
157  int p1; _ArrayIndex1DWO_(solver->ndims,dim,indexb,offset,ghosts,p1);
158  int p2; _ArrayIndex1D_ (solver->ndims,dim,indexi,ghosts,p2);
159  _ArrayCopy1D_((adv+p2*adv_nvar), (adv+p1*adv_nvar), adv_nvar);
160  _ArrayIncrementIndex_(solver->ndims,bounds,indexb,done);
161  }
162  } else {
163  /* left boundary */
164  if (!mpi->ip[d]) {
165  _ArrayCopy1D_(dim,bounds,solver->ndims); bounds[d] = ghosts;
166  _ArraySetValue_(offset,solver->ndims,0); offset[d] = -ghosts;
167  done = 0; _ArraySetValue_(indexb,solver->ndims,0);
168  while (!done) {
169  _ArrayCopy1D_(indexb,indexi,solver->ndims); indexi[d] = ghosts-1-indexb[d];
170  int p1; _ArrayIndex1DWO_(solver->ndims,dim,indexb,offset,ghosts,p1);
171  int p2; _ArrayIndex1D_ (solver->ndims,dim,indexi,ghosts,p2);
172  _ArrayCopy1D_((adv+p2*adv_nvar), (adv+p1*adv_nvar), adv_nvar);
173  _ArrayIncrementIndex_(solver->ndims,bounds,indexb,done);
174  }
175  }
176  /* right boundary */
177  if (mpi->ip[d] == mpi->iproc[d]-1) {
178  _ArrayCopy1D_(dim,bounds,solver->ndims); bounds[d] = ghosts;
179  _ArraySetValue_(offset,solver->ndims,0); offset[d] = dim[d];
180  done = 0; _ArraySetValue_(indexb,solver->ndims,0);
181  while (!done) {
182  _ArrayCopy1D_(indexb,indexi,solver->ndims); indexi[d] = dim[d]-1-indexb[d];
183  int p1; _ArrayIndex1DWO_(solver->ndims,dim,indexb,offset,ghosts,p1);
184  int p2; _ArrayIndex1D_ (solver->ndims,dim,indexi,ghosts,p2);
185  _ArrayCopy1D_((adv+p2*adv_nvar), (adv+p1*adv_nvar), adv_nvar);
186  _ArrayIncrementIndex_(solver->ndims,bounds,indexb,done);
187  }
188  }
189  }
190  }
191 
192  return(0);
193 }
int npoints_local_wghosts
Definition: hypar.h:42
#define _ArraySetValue_(x, size, value)
char adv_filename[_MAX_STRING_SIZE_]
Definition: linearadr.h:43
int * dim_global
Definition: hypar.h:33
#define _ArrayIncrementIndex_(N, imax, i, done)
void GetStringFromInteger(int, char *, int)
int adv_arr_size
Definition: linearadr.h:47
void * physics
Definition: hypar.h:266
Structure containing variables and parameters specific to the linear advection-diffusion-reaction mod...
Definition: linearadr.h:37
int * dim_local
Definition: hypar.h:37
double * a
Definition: linearadr.h:50
int ghosts
Definition: hypar.h:52
#define _ArrayIndex1D_(N, imax, i, ghost, index)
int ReadArray(int, int, int *, int *, int, void *, void *, double *, double *, char *, int *)
Definition: ReadArray.c:25
#define _MAX_STRING_SIZE_
Definition: basic.h:14
int * isPeriodic
Definition: hypar.h:162
#define _ArrayIndex1DWO_(N, imax, i, offset, ghost, index)
#define _ArrayCopy1D_(x, y, size)
int nvars
Definition: hypar.h:29
int ReadArraywInterp(int, int, int *, int *, int *, int, void *, void *, double *, double *, char *, int *)
int constant_advection
Definition: linearadr.h:40
int ndims
Definition: hypar.h:26
#define IERR
Definition: basic.h:16
int MPIExchangeBoundariesnD(int, int, int *, int, void *, double *)
Structure of MPI-related variables.
Structure containing all solver-specific variables and functions.
Definition: hypar.h:23