Commit 80229ccb authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

added scrip_remap_conserv.c

parent caf9c668
......@@ -591,10 +591,13 @@ src/realtime.c -text
src/remap.h -text
src/remap_scrip_io.c -text
src/remap_search_latbins.c -text
src/remap_store_link.c -text
src/remap_store_link.h -text
src/remaplib.c -text
src/remapsort.c -text
src/results_template_parser.c -text
src/results_template_parser.h -text
src/scrip_remap_conserv.c -text
src/specspace.c -text
src/specspace.h -text
src/statistic.c -text
......
......@@ -271,6 +271,9 @@ cdo_SOURCES += Adisit.c \
remapsort.c \
remap_scrip_io.c \
remap_search_latbins.c \
remap_store_link.c \
remap_store_link.h \
scrip_remap_conserv.c \
stdnametable.c \
stdnametable.h \
specspace.c \
......
......@@ -163,6 +163,7 @@ am__cdo_SOURCES_DIST = cdo.c Adisit.c Arith.c Arithc.c Arithdays.c \
pstream.c pstream.h pstream_int.h pthread_debug.c \
pthread_debug.h readline.c realtime.c remap.h remaplib.c \
remapsort.c remap_scrip_io.c remap_search_latbins.c \
remap_store_link.c remap_store_link.h scrip_remap_conserv.c \
stdnametable.c stdnametable.h specspace.c specspace.h \
statistic.c statistic.h table.c text.c text.h timebase.h \
timer.c userlog.c util.c util.h vinterp.c vinterp.h zaxis.c \
......@@ -292,7 +293,9 @@ am_cdo_OBJECTS = cdo-cdo.$(OBJEXT) cdo-Adisit.$(OBJEXT) \
cdo-pthread_debug.$(OBJEXT) cdo-readline.$(OBJEXT) \
cdo-realtime.$(OBJEXT) cdo-remaplib.$(OBJEXT) \
cdo-remapsort.$(OBJEXT) cdo-remap_scrip_io.$(OBJEXT) \
cdo-remap_search_latbins.$(OBJEXT) cdo-stdnametable.$(OBJEXT) \
cdo-remap_search_latbins.$(OBJEXT) \
cdo-remap_store_link.$(OBJEXT) \
cdo-scrip_remap_conserv.$(OBJEXT) cdo-stdnametable.$(OBJEXT) \
cdo-specspace.$(OBJEXT) cdo-statistic.$(OBJEXT) \
cdo-table.$(OBJEXT) cdo-text.$(OBJEXT) cdo-timer.$(OBJEXT) \
cdo-userlog.$(OBJEXT) cdo-util.$(OBJEXT) cdo-vinterp.$(OBJEXT) \
......@@ -594,7 +597,8 @@ cdo_SOURCES = cdo.c Adisit.c Arith.c Arithc.c Arithdays.c Arithlat.c \
pipe.c pipe.h printinfo.h process.c process.h pstream.c \
pstream.h pstream_int.h pthread_debug.c pthread_debug.h \
readline.c realtime.c remap.h remaplib.c remapsort.c \
remap_scrip_io.c remap_search_latbins.c stdnametable.c \
remap_scrip_io.c remap_search_latbins.c remap_store_link.c \
remap_store_link.h scrip_remap_conserv.c stdnametable.c \
stdnametable.h specspace.c specspace.h statistic.c statistic.h \
table.c text.c text.h timebase.h timer.c userlog.c util.c \
util.h vinterp.c vinterp.h zaxis.c clipping/clipping.c \
......@@ -984,9 +988,11 @@ distclean-compile:
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-realtime.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-remap_scrip_io.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-remap_search_latbins.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-remap_store_link.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-remaplib.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-remapsort.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-results_template_parser.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-scrip_remap_conserv.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-specspace.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-statistic.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/cdo-stdnametable.Po@am__quote@
......@@ -4256,6 +4262,34 @@ cdo-remap_search_latbins.obj: remap_search_latbins.c
@AMDEP_TRUE@@am__fastdepCC_FALSE@ DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
@am__fastdepCC_FALSE@ $(AM_V_CC@am__nodep@)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -c -o cdo-remap_search_latbins.obj `if test -f 'remap_search_latbins.c'; then $(CYGPATH_W) 'remap_search_latbins.c'; else $(CYGPATH_W) '$(srcdir)/remap_search_latbins.c'; fi`
cdo-remap_store_link.o: remap_store_link.c
@am__fastdepCC_TRUE@ $(AM_V_CC)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -MT cdo-remap_store_link.o -MD -MP -MF $(DEPDIR)/cdo-remap_store_link.Tpo -c -o cdo-remap_store_link.o `test -f 'remap_store_link.c' || echo '$(srcdir)/'`remap_store_link.c
@am__fastdepCC_TRUE@ $(AM_V_at)$(am__mv) $(DEPDIR)/cdo-remap_store_link.Tpo $(DEPDIR)/cdo-remap_store_link.Po
@AMDEP_TRUE@@am__fastdepCC_FALSE@ $(AM_V_CC)source='remap_store_link.c' object='cdo-remap_store_link.o' libtool=no @AMDEPBACKSLASH@
@AMDEP_TRUE@@am__fastdepCC_FALSE@ DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
@am__fastdepCC_FALSE@ $(AM_V_CC@am__nodep@)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -c -o cdo-remap_store_link.o `test -f 'remap_store_link.c' || echo '$(srcdir)/'`remap_store_link.c
cdo-remap_store_link.obj: remap_store_link.c
@am__fastdepCC_TRUE@ $(AM_V_CC)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -MT cdo-remap_store_link.obj -MD -MP -MF $(DEPDIR)/cdo-remap_store_link.Tpo -c -o cdo-remap_store_link.obj `if test -f 'remap_store_link.c'; then $(CYGPATH_W) 'remap_store_link.c'; else $(CYGPATH_W) '$(srcdir)/remap_store_link.c'; fi`
@am__fastdepCC_TRUE@ $(AM_V_at)$(am__mv) $(DEPDIR)/cdo-remap_store_link.Tpo $(DEPDIR)/cdo-remap_store_link.Po
@AMDEP_TRUE@@am__fastdepCC_FALSE@ $(AM_V_CC)source='remap_store_link.c' object='cdo-remap_store_link.obj' libtool=no @AMDEPBACKSLASH@
@AMDEP_TRUE@@am__fastdepCC_FALSE@ DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
@am__fastdepCC_FALSE@ $(AM_V_CC@am__nodep@)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -c -o cdo-remap_store_link.obj `if test -f 'remap_store_link.c'; then $(CYGPATH_W) 'remap_store_link.c'; else $(CYGPATH_W) '$(srcdir)/remap_store_link.c'; fi`
cdo-scrip_remap_conserv.o: scrip_remap_conserv.c
@am__fastdepCC_TRUE@ $(AM_V_CC)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -MT cdo-scrip_remap_conserv.o -MD -MP -MF $(DEPDIR)/cdo-scrip_remap_conserv.Tpo -c -o cdo-scrip_remap_conserv.o `test -f 'scrip_remap_conserv.c' || echo '$(srcdir)/'`scrip_remap_conserv.c
@am__fastdepCC_TRUE@ $(AM_V_at)$(am__mv) $(DEPDIR)/cdo-scrip_remap_conserv.Tpo $(DEPDIR)/cdo-scrip_remap_conserv.Po
@AMDEP_TRUE@@am__fastdepCC_FALSE@ $(AM_V_CC)source='scrip_remap_conserv.c' object='cdo-scrip_remap_conserv.o' libtool=no @AMDEPBACKSLASH@
@AMDEP_TRUE@@am__fastdepCC_FALSE@ DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
@am__fastdepCC_FALSE@ $(AM_V_CC@am__nodep@)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -c -o cdo-scrip_remap_conserv.o `test -f 'scrip_remap_conserv.c' || echo '$(srcdir)/'`scrip_remap_conserv.c
cdo-scrip_remap_conserv.obj: scrip_remap_conserv.c
@am__fastdepCC_TRUE@ $(AM_V_CC)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -MT cdo-scrip_remap_conserv.obj -MD -MP -MF $(DEPDIR)/cdo-scrip_remap_conserv.Tpo -c -o cdo-scrip_remap_conserv.obj `if test -f 'scrip_remap_conserv.c'; then $(CYGPATH_W) 'scrip_remap_conserv.c'; else $(CYGPATH_W) '$(srcdir)/scrip_remap_conserv.c'; fi`
@am__fastdepCC_TRUE@ $(AM_V_at)$(am__mv) $(DEPDIR)/cdo-scrip_remap_conserv.Tpo $(DEPDIR)/cdo-scrip_remap_conserv.Po
@AMDEP_TRUE@@am__fastdepCC_FALSE@ $(AM_V_CC)source='scrip_remap_conserv.c' object='cdo-scrip_remap_conserv.obj' libtool=no @AMDEPBACKSLASH@
@AMDEP_TRUE@@am__fastdepCC_FALSE@ DEPDIR=$(DEPDIR) $(CCDEPMODE) $(depcomp) @AMDEPBACKSLASH@
@am__fastdepCC_FALSE@ $(AM_V_CC@am__nodep@)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -c -o cdo-scrip_remap_conserv.obj `if test -f 'scrip_remap_conserv.c'; then $(CYGPATH_W) 'scrip_remap_conserv.c'; else $(CYGPATH_W) '$(srcdir)/scrip_remap_conserv.c'; fi`
cdo-stdnametable.o: stdnametable.c
@am__fastdepCC_TRUE@ $(AM_V_CC)$(CC) $(DEFS) $(DEFAULT_INCLUDES) $(INCLUDES) $(cdo_CPPFLAGS) $(CPPFLAGS) $(AM_CFLAGS) $(CFLAGS) -MT cdo-stdnametable.o -MD -MP -MF $(DEPDIR)/cdo-stdnametable.Tpo -c -o cdo-stdnametable.o `test -f 'stdnametable.c' || echo '$(srcdir)/'`stdnametable.c
@am__fastdepCC_TRUE@ $(AM_V_at)$(am__mv) $(DEPDIR)/cdo-stdnametable.Tpo $(DEPDIR)/cdo-stdnametable.Po
......
......@@ -208,7 +208,6 @@ void print_remap_info(int operfunc, int map_type, remapgrid_t *src_grid, remapgr
double remap_threshhold = 2;
int remap_test = 0;
int remap_order = 1;
int remap_store_link_fast = TRUE;
int remap_non_global = FALSE;
int remap_num_srch_bins = 180;
int lremap_num_srch_bins = FALSE;
......@@ -356,6 +355,8 @@ void get_remap_env(void)
}
}
int remap_store_link_fast = TRUE;
envstr = getenv("REMAP_STORE_LINK_FAST");
if ( envstr )
{
......
#ifndef _REMAP_H
#define _REMAP_H
#if defined(_OPENMP)
#include <omp.h>
#endif
......@@ -8,6 +11,20 @@
#define PI2 2.0*PI
#define PIH 0.5*PI
#define ZERO 0.0
#define ONE 1.0
#define TWO 2.0
#define THREE 3.0
#define HALF 0.5
#define QUART 0.25
#define BIGNUM 1.e+20
#define TINY 1.e-14
#define REMAP_GRID_TYPE_REG2D 1
#define REMAP_GRID_TYPE_CURVE2D 2
#define REMAP_GRID_TYPE_UNSTRUCT 3
#define REMAP_GRID_BASIS_SRC 1
#define REMAP_GRID_BASIS_TGT 2
......@@ -194,3 +211,8 @@ void store_link_bilin(remapvars_t *rv, int dst_add, int *restrict src_add, doubl
void calc_bin_addr(long gridsize, long nbins, const restr_t* restrict bin_lats, const restr_t* restrict cell_bound_box, int* restrict bin_addr);
void calc_lat_bins(remapgrid_t* src_grid, remapgrid_t* tgt_grid, int map_type);
long get_srch_cells(long tgt_grid_add, long nbins, int *bin_addr1, int *bin_addr2,
restr_t *tgt_cell_bound_box, restr_t *src_cell_bound_box, long src_grid_size, int *srch_add);
#endif /* _REMAP_H */
#include "cdo.h"
#include "remap.h"
void calc_bin_addr(long gridsize, long nbins, const restr_t* restrict bin_lats, const restr_t* restrict cell_bound_box, int* restrict bin_addr)
{
long n, n2, nele, nele4;
......@@ -123,3 +124,92 @@ void calc_lat_bins(remapgrid_t* src_grid, remapgrid_t* tgt_grid, int map_type)
free(src_grid->cell_bound_box); src_grid->cell_bound_box = NULL;
}
}
long get_srch_cells(long tgt_grid_add, long nbins, int *bin_addr1, int *bin_addr2,
restr_t *tgt_cell_bound_box, restr_t *src_cell_bound_box, long src_grid_size, int *srch_add)
{
long num_srch_cells; /* num cells in restricted search arrays */
long min_add; /* addresses for restricting search of */
long max_add; /* destination grid */
long n, n2; /* generic counters */
long src_grid_add; /* current linear address for src cell */
long src_grid_addm4;
restr_t bound_box_lat1, bound_box_lat2, bound_box_lon1, bound_box_lon2;
/* Restrict searches first using search bins */
min_add = src_grid_size - 1;
max_add = 0;
for ( n = 0; n < nbins; ++n )
{
n2 = n<<1;
if ( tgt_grid_add >= bin_addr1[n2] && tgt_grid_add <= bin_addr1[n2+1] )
{
if ( bin_addr2[n2 ] < min_add ) min_add = bin_addr2[n2 ];
if ( bin_addr2[n2+1] > max_add ) max_add = bin_addr2[n2+1];
}
}
/* Further restrict searches using bounding boxes */
bound_box_lat1 = tgt_cell_bound_box[0];
bound_box_lat2 = tgt_cell_bound_box[1];
bound_box_lon1 = tgt_cell_bound_box[2];
bound_box_lon2 = tgt_cell_bound_box[3];
num_srch_cells = 0;
for ( src_grid_add = min_add; src_grid_add <= max_add; ++src_grid_add )
{
src_grid_addm4 = src_grid_add<<2;
if ( (src_cell_bound_box[src_grid_addm4+2] <= bound_box_lon2) &&
(src_cell_bound_box[src_grid_addm4+3] >= bound_box_lon1) )
{
if ( (src_cell_bound_box[src_grid_addm4 ] <= bound_box_lat2) &&
(src_cell_bound_box[src_grid_addm4+1] >= bound_box_lat1) )
{
srch_add[num_srch_cells] = src_grid_add;
num_srch_cells++;
}
}
}
if ( bound_box_lon1 < RESTR_SCALE(0.) || bound_box_lon2 > RESTR_SCALE(PI2) )
{
if ( bound_box_lon1 < RESTR_SCALE(0.) )
{
bound_box_lon1 += RESTR_SCALE(PI2);
bound_box_lon2 += RESTR_SCALE(PI2);
}
else
{
bound_box_lon1 -= RESTR_SCALE(PI2);
bound_box_lon2 -= RESTR_SCALE(PI2);
}
for ( src_grid_add = min_add; src_grid_add <= max_add; ++src_grid_add )
{
src_grid_addm4 = src_grid_add<<2;
if ( (src_cell_bound_box[src_grid_addm4+2] <= bound_box_lon2) &&
(src_cell_bound_box[src_grid_addm4+3] >= bound_box_lon1) )
{
if ( (src_cell_bound_box[src_grid_addm4 ] <= bound_box_lat2) &&
(src_cell_bound_box[src_grid_addm4+1] >= bound_box_lat1) )
{
long ii;
for ( ii = 0; ii < num_srch_cells; ++ii )
if ( srch_add[ii] == src_grid_add ) break;
if ( ii == num_srch_cells )
{
srch_add[num_srch_cells] = src_grid_add;
num_srch_cells++;
}
}
}
}
}
return (num_srch_cells);
}
#include "cdo.h"
#include "cdo_int.h"
#include "remap.h"
#include "remap_store_link.h"
int remap_store_link_fast = TRUE;
void grid_store_init(grid_store_t* grid_store, long gridsize)
{
long iblk;
long blksize[] = {128, 256, 512, 1024, 2048, 4096, 8192};
long nblks = sizeof(blksize)/sizeof(long);
long nblocks;
for ( iblk = nblks-1; iblk >= 0; --iblk )
if ( gridsize/blksize[iblk] > 99 ) break;
if ( iblk < 0 ) iblk = 0;
/* grid_store->blk_size = BLK_SIZE; */
grid_store->blk_size = blksize[iblk];
grid_store->max_size = gridsize;
grid_store->nblocks = grid_store->max_size/grid_store->blk_size;
if ( grid_store->max_size%grid_store->blk_size > 0 ) grid_store->nblocks++;
if ( cdoVerbose )
fprintf(stdout, "blksize = %d lastblksize = %d max_size = %d nblocks = %d\n",
grid_store->blk_size, grid_store->max_size%grid_store->blk_size,
grid_store->max_size, grid_store->nblocks);
grid_store->blksize = malloc(grid_store->nblocks*sizeof(int));
grid_store->nlayers = malloc(grid_store->nblocks*sizeof(int));
grid_store->layers = malloc(grid_store->nblocks*sizeof(grid_layer_t *));
nblocks = grid_store->nblocks;
for ( iblk = 0; iblk < nblocks; ++iblk )
{
grid_store->blksize[iblk] = grid_store->blk_size;
grid_store->nlayers[iblk] = 0;
grid_store->layers[iblk] = NULL;
}
if ( grid_store->max_size%grid_store->blk_size > 0 )
grid_store->blksize[grid_store->nblocks-1] = grid_store->max_size%grid_store->blk_size;
}
void grid_store_delete(grid_store_t* grid_store)
{
grid_layer_t *grid_layer, *grid_layer_f;
long ilayer;
long i, j;
long iblk;
for ( iblk = 0; iblk < grid_store->nblocks; ++iblk )
{
j = 0;
grid_layer = grid_store->layers[iblk];
for ( ilayer = 0; ilayer < grid_store->nlayers[iblk]; ++ilayer )
{
if ( cdoVerbose )
{
for ( i = 0; i < grid_store->blksize[iblk]; ++i )
if ( grid_layer->grid2_link[i] != -1 ) j++;
}
grid_layer_f = grid_layer;
free(grid_layer->grid2_link);
grid_layer = grid_layer->next;
free(grid_layer_f);
}
if ( cdoVerbose )
{
fprintf(stderr, "block = %ld nlayers = %d allocated = %d used = %ld\n",
iblk+1, grid_store->nlayers[iblk],
grid_store->nlayers[iblk]*grid_store->blksize[iblk], j);
}
}
free(grid_store->blksize);
free(grid_store->layers);
free(grid_store->nlayers);
}
/*
This routine stores the address and weight for this link in the appropriate
address and weight arrays and resizes those arrays if necessary.
*/
void store_link_cnsrv_fast(remapvars_t* rv, long add1, long add2, long num_wts, double* weights, grid_store_t* grid_store)
{
/*
Input variables:
int add1 ! address on source grid
int add2 ! address on target grid
double weights[] ! array of remapping weights for this link
*/
/* Local variables */
long nlink; /* link index */
long ilayer, i, iblk, iadd2;
long nlayer, blksize;
int lstore_link;
grid_layer_t *grid_layer, **grid_layer2;
/* If all weights are ZERO, do not bother storing the link */
if ( num_wts == 3 )
{
if ( IS_EQUAL(weights[0], 0) && IS_EQUAL(weights[1], 0) && IS_EQUAL(weights[2], 0) ) return;
}
else
{
if ( IS_EQUAL(weights[0], 0) ) return;
}
/* If the link already exists, add the weight to the current weight arrays */
iblk = BLK_NUM(add2);
iadd2 = BLK_IDX(add2);
lstore_link = FALSE;
grid_layer2 = &grid_store->layers[iblk];
nlayer = grid_store->nlayers[iblk];
for ( ilayer = 0; ilayer < nlayer; ++ilayer )
{
grid_layer = *grid_layer2;
nlink = grid_layer->grid2_link[iadd2];
if ( nlink == -1 )
{
break;
}
else if ( add1 == rv->src_grid_add[nlink] )
{
lstore_link = TRUE;
break;
}
grid_layer2 = &(*grid_layer2)->next;
}
if ( lstore_link )
{
for ( i = 0; i < num_wts; ++i ) rv->wts[num_wts*nlink+i] += weights[i];
return;
}
/*
If the link does not yet exist, increment number of links and
check to see if remap arrays need to be increased to accomodate
the new link. Then store the link.
*/
nlink = rv->num_links;
if ( ilayer < grid_store->nlayers[iblk] )
{
grid_layer->grid2_link[iadd2] = nlink;
}
else
{
grid_layer = malloc(sizeof(grid_layer_t));
grid_layer->next = NULL;
grid_layer->grid2_link = malloc(grid_store->blksize[iblk]*sizeof(int));
blksize = grid_store->blksize[iblk];
for ( i = 0; i < blksize; ++i )
grid_layer->grid2_link[i] = -1;
grid_layer->grid2_link[iadd2] = nlink;
*grid_layer2 = grid_layer;
grid_store->nlayers[iblk]++;
}
rv->num_links++;
if ( rv->num_links >= rv->max_links )
resize_remap_vars(rv, rv->resize_increment);
rv->src_grid_add[nlink] = add1;
rv->tgt_grid_add[nlink] = add2;
for ( i = 0; i < num_wts; ++i ) rv->wts[num_wts*nlink+i] = weights[i];
} /* store_link_cnsrv_fast */
/*
This routine stores the address and weight for this link in the appropriate
address and weight arrays and resizes those arrays if necessary.
*/
void store_link_cnsrv(remapvars_t *rv, long add1, long add2, double *restrict weights, int *link_add1[2], int *link_add2[2])
{
/*
Input variables:
int add1 ! address on source grid
int add2 ! address on target grid
double weights[3] ! array of remapping weights for this link
*/
/* Local variables */
long nlink, min_link, max_link; /* link index */
/* If all weights are ZERO, do not bother storing the link */
if ( IS_EQUAL(weights[0], 0) && IS_EQUAL(weights[1], 0) && IS_EQUAL(weights[2], 0) ) return;
/* Restrict the range of links to search for existing links */
min_link = MIN(link_add1[0][add1], link_add2[0][add2]);
max_link = MAX(link_add1[1][add1], link_add2[1][add2]);
if ( min_link == -1 )
{
min_link = 0;
max_link = -1;
}
/* If the link already exists, add the weight to the current weight arrays */
#if defined(SX)
#define STRIPED 1
#if STRIPED
#define STRIPLENGTH 4096
{
long ilink = max_link + 1;
long strip, estrip;
nlink = 0;
for ( strip=min_link; strip <= max_link; strip+=STRIPLENGTH )
{
estrip = MIN(max_link-strip+1, STRIPLENGTH);
for ( nlink = 0; nlink < estrip; ++nlink )
{
if ( add2 == rv->tgt_grid_add[strip+nlink] &&
add1 == rv->src_grid_add[strip+nlink] )
ilink = strip + nlink;
}
if (ilink != (max_link + 1)) break;
}
nlink += strip;
if (ilink != (max_link + 1)) nlink = ilink;
}
#else
{
long ilink = max_link + 1;
for ( nlink = min_link; nlink <= max_link; ++nlink )
{
if ( add2 == rv->tgt_grid_add[nlink] )
if ( add1 == rv->src_grid_add[nlink] ) ilink = nlink;
}
if ( ilink != (max_link + 1) ) nlink = ilink;
}
#endif
#else
for ( nlink = min_link; nlink <= max_link; ++nlink )
{
if ( add2 == rv->tgt_grid_add[nlink] )
if ( add1 == rv->src_grid_add[nlink] ) break;
}
#endif
if ( nlink <= max_link )
{
rv->wts[3*nlink ] += weights[0];
rv->wts[3*nlink+1] += weights[1];
rv->wts[3*nlink+2] += weights[2];
return;
}
/*
If the link does not yet exist, increment number of links and
check to see if remap arrays need to be increased to accomodate
the new link. Then store the link.
*/
nlink = rv->num_links;
rv->num_links++;
if ( rv->num_links >= rv->max_links )
resize_remap_vars(rv, rv->resize_increment);
rv->src_grid_add[nlink] = add1;
rv->tgt_grid_add[nlink] = add2;
rv->wts[3*nlink ] = weights[0];
rv->wts[3*nlink+1] = weights[1];
rv->wts[3*nlink+2] = weights[2];
if ( link_add1[0][add1] == -1 ) link_add1[0][add1] = (int)nlink;
if ( link_add2[0][add2] == -1 ) link_add2[0][add2] = (int)nlink;
link_add1[1][add1] = (int)nlink;
link_add2[1][add2] = (int)nlink;
} /* store_link_cnsrv */
#ifndef _REMAP_STORE_LINK_H
#define _REMAP_STORE_LINK_H
extern int remap_store_link_fast;
/* used for store_link_fast */
#define BLK_SIZE 4096
#define BLK_NUM(x) (x/grid_store->blk_size)
#define BLK_IDX(x) (x%grid_store->blk_size)
struct grid_layer
{
int *grid2_link;
struct grid_layer *next;
};
typedef struct grid_layer grid_layer_t;
typedef struct
{
int blk_size;
int max_size;
int nblocks;
int *blksize;
int *nlayers;
grid_layer_t **layers;
} grid_store_t;
void grid_store_init(grid_store_t *grid_store, long gridsize);
void grid_store_delete(grid_store_t *grid_store);
void store_link_cnsrv_fast(remapvars_t *rv, long add1, long add2, long num_wts, double *weights, grid_store_t *grid_store);
void store_link_cnsrv(remapvars_t *rv, long add1, long add2, double *restrict weights, int *link_add1[2], int *link_add2[2]);
#endif /* _REMAP_STORE_LINK */
This diff is collapsed.
This diff is collapsed.
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment