RTOp_TOp_set_sub_vector.c

00001 /*
00002 // @HEADER
00003 // ***********************************************************************
00004 // 
00005 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
00006 //                  Copyright (2003) Sandia Corporation
00007 // 
00008 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 // license for use of this work by or on behalf of the U.S. Government.
00010 // 
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //  
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //  
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA
00025 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
00026 // 
00027 // ***********************************************************************
00028 // @HEADER
00029 */
00030 
00031 #include "RTOp_TOp_set_sub_vector.h"
00032 #include "RTOp_obj_null_vtbl.h"
00033 
00034 /* Operator object data virtual function table */
00035 
00036 struct RTOp_TOp_set_sub_vector_state_t { /* operator object instance data */
00037   struct RTOp_SparseSubVector  sub_vec;   /* The sub vector to set! */
00038   int                    owns_mem;  /* if true then memory is sub_vec must be deleted! */
00039 };
00040 
00041 static int get_op_type_num_entries(
00042   const struct RTOp_obj_type_vtbl_t*   vtbl
00043   ,const void* obj_data
00044   ,int* num_values
00045   ,int* num_indexes
00046   ,int* num_chars
00047   )
00048 {
00049   const struct RTOp_TOp_set_sub_vector_state_t *state = NULL;
00050 #ifdef RTOp_DEBUG
00051   assert( num_values );
00052   assert( num_indexes );
00053   assert( num_chars );
00054   assert( obj_data );
00055 #endif
00056   state = (const struct RTOp_TOp_set_sub_vector_state_t*)obj_data;
00057   *num_values  = state->sub_vec.sub_nz;  /* values[] */
00058   *num_indexes =
00059     8  /* global_offset, sub_dim, sub_nz, values_stride, indices_stride, local_offset, is_sorted, owns_mem */
00060     + (state->sub_vec.indices ? state->sub_vec.sub_nz : 0 ); /* indices[] */
00061   *num_chars   = 0;
00062   return 0;
00063 }
00064 
00065 static int op_create(
00066   const struct RTOp_obj_type_vtbl_t* vtbl, const void* instance_data
00067   , RTOp_ReductTarget* obj )
00068 {
00069   struct RTOp_TOp_set_sub_vector_state_t *state = NULL;
00070   *obj = malloc(sizeof(struct RTOp_TOp_set_sub_vector_state_t));
00071   state = (struct RTOp_TOp_set_sub_vector_state_t*)*obj;
00072   RTOp_sparse_sub_vector_null( &state->sub_vec );
00073   state->owns_mem = 0;
00074   return 0;
00075 }
00076 
00077 static int op_free(
00078   const struct RTOp_obj_type_vtbl_t* vtbl, const void* dummy
00079   , void ** obj_data )
00080 {
00081   struct RTOp_TOp_set_sub_vector_state_t *state = NULL;
00082 #ifdef RTOp_DEBUG
00083   assert( *obj_data );
00084 #endif
00085   state = (struct RTOp_TOp_set_sub_vector_state_t*)*obj_data;
00086   if( state->owns_mem ) {
00087     free( (void*)state->sub_vec.values );
00088     free( (void*)state->sub_vec.indices );
00089   }
00090   free( *obj_data );
00091   *obj_data = NULL;
00092   return 0;
00093 }
00094 
00095 static int extract_op_state(
00096   const struct RTOp_obj_type_vtbl_t*   vtbl
00097   ,const void *       dummy
00098   ,void *             obj_data
00099   ,int                num_values
00100   ,RTOp_value_type    value_data[]
00101   ,int                num_indexes
00102   ,RTOp_index_type    index_data[]
00103   ,int                num_chars
00104   ,RTOp_char_type     char_data[]
00105   )
00106 {
00107   const struct RTOp_TOp_set_sub_vector_state_t *state = NULL;
00108   register RTOp_index_type k, j;
00109 #ifdef RTOp_DEBUG
00110   assert( obj_data );
00111 #endif
00112   state = (const struct RTOp_TOp_set_sub_vector_state_t*)obj_data;
00113 #ifdef RTOp_DEBUG
00114   assert( num_values == state->sub_vec.sub_nz );
00115   assert( num_indexes == 8 + (state->sub_vec.indices ? state->sub_vec.sub_nz : 0 ) );
00116   assert( num_chars == 0 );
00117   assert( value_data );
00118   assert( index_data );
00119 #endif
00120   for( k = 0; k < state->sub_vec.sub_nz; ++k )
00121     value_data[k] = state->sub_vec.values[k];
00122   index_data[k=0] = state->sub_vec.global_offset;
00123   index_data[++k] = state->sub_vec.sub_dim;
00124   index_data[++k] = state->sub_vec.sub_nz;
00125   index_data[++k] = state->sub_vec.values_stride;
00126   index_data[++k] = state->sub_vec.indices_stride;
00127   index_data[++k] = state->sub_vec.local_offset;
00128   index_data[++k] = state->sub_vec.is_sorted;
00129   index_data[++k] = state->owns_mem;
00130   if( state->sub_vec.indices ) {
00131     for( j = 0; j < state->sub_vec.sub_nz; ++j )
00132       index_data[++k] = state->sub_vec.indices[j];
00133   }
00134   return 0;
00135 }
00136 
00137 static int load_op_state(
00138   const struct RTOp_obj_type_vtbl_t*   vtbl
00139   ,const void *            dummy
00140   ,int                     num_values
00141   ,const RTOp_value_type   value_data[]
00142   ,int                     num_indexes
00143   ,const RTOp_index_type   index_data[]
00144   ,int                     num_chars
00145   ,const RTOp_char_type    char_data[]
00146   ,void **                 obj_data
00147   )
00148 {
00149   struct RTOp_TOp_set_sub_vector_state_t *state = NULL;
00150   register RTOp_index_type k, j;
00151   RTOp_value_type *values = NULL;
00152   /* Allocate the operator object's state data if it has not been already */
00153   if( *obj_data == NULL ) {
00154     *obj_data = (void*)malloc(sizeof(struct RTOp_TOp_set_sub_vector_state_t));
00155     state = (struct RTOp_TOp_set_sub_vector_state_t*)*obj_data;
00156     RTOp_sparse_sub_vector_null( &state->sub_vec );
00157     state->owns_mem = 0;
00158   }
00159   /* Get the operator object's state data */
00160 #ifdef RTOp_DEBUG
00161   assert( *obj_data );
00162 #endif
00163   state = (struct RTOp_TOp_set_sub_vector_state_t*)*obj_data;
00164 #ifdef RTOp_DEBUG
00165   if( num_indexes > 8 ) assert( num_values == num_indexes - 8 );
00166 #endif
00167   /* Reallocate storage if we have to */
00168   if( num_values != state->sub_vec.sub_nz || !state->owns_mem ) {
00169     /* Delete current storage if owned */
00170     if( state->owns_mem ) {
00171       free( (RTOp_value_type*)state->sub_vec.values );
00172       free( (RTOp_index_type*)state->sub_vec.indices );
00173     }
00174     /* We need to reallocate storage for values[] and perhaps */
00175     state->sub_vec.values = (RTOp_value_type*)malloc(num_values*sizeof(RTOp_value_type));
00176     if( num_indexes > 8 )
00177       state->sub_vec.indices = (RTOp_index_type*)malloc(num_values*sizeof(RTOp_index_type));
00178     state->owns_mem = 1;
00179   }
00180   /* Set the internal state! */
00181 #ifdef RTOp_DEBUG
00182   assert( num_chars == 0 );
00183   assert( value_data );
00184   assert( index_data );
00185 #endif
00186   for( values = (RTOp_value_type*)state->sub_vec.values, k = 0; k < num_values; ++k )
00187     *values++ = value_data[k];
00188   state->sub_vec.global_offset  = index_data[k=0];
00189   state->sub_vec.sub_dim        = index_data[++k];
00190   state->sub_vec.sub_nz         = index_data[++k];
00191   state->sub_vec.values_stride  = index_data[++k];
00192   state->sub_vec.indices_stride = index_data[++k];
00193   state->sub_vec.local_offset   = index_data[++k];
00194   state->sub_vec.is_sorted      = index_data[++k];
00195   state->owns_mem               = index_data[++k];
00196   if( num_indexes > 8 ) {
00197     for( j = 0; j < num_values; ++j )
00198       ((RTOp_index_type*)state->sub_vec.indices)[j] = index_data[++k];
00199   }
00200   return 0;
00201 }
00202 
00203 static struct RTOp_obj_type_vtbl_t  instance_obj_vtbl =
00204 {
00205   get_op_type_num_entries
00206   ,op_create
00207   ,NULL
00208   ,op_free
00209   ,extract_op_state
00210   ,load_op_state
00211 };
00212 
00213 /* Implementation functions */
00214 
00215 static int RTOp_TOp_set_sub_vector_apply_op(
00216   const struct RTOp_RTOp_vtbl_t* vtbl, const void* obj_data
00217   , const int num_vecs, const struct RTOp_SubVector vecs[]
00218   , const int num_targ_vecs, const struct RTOp_MutableSubVector targ_vecs[]
00219   , RTOp_ReductTarget targ_obj )
00220 {
00221   const struct RTOp_TOp_set_sub_vector_state_t *state = NULL;
00222   RTOp_index_type        v_global_offset;
00223   RTOp_index_type        v_sub_dim;
00224   RTOp_index_type        v_sub_nz;
00225   const RTOp_value_type  *v_val;
00226   ptrdiff_t              v_val_s;
00227   const RTOp_index_type  *v_ind;
00228   ptrdiff_t              v_ind_s;
00229   ptrdiff_t              v_l_off;
00230   int                    v_sorted;
00231   RTOp_index_type        z_global_offset;
00232   RTOp_index_type        z_sub_dim;
00233   RTOp_value_type        *z_val;
00234   ptrdiff_t              z_val_s;
00235   register RTOp_index_type  k, i;
00236   RTOp_index_type num_overlap;
00237 
00238   /* */
00239   /* Validate the input */
00240   /* */
00241   if( num_vecs != 0 )
00242     return RTOp_ERR_INVALID_NUM_VECS;
00243   if( num_targ_vecs != 1 )
00244     return RTOp_ERR_INVALID_NUM_TARG_VECS;
00245   assert(targ_vecs);
00246   assert(targ_obj == RTOp_REDUCT_OBJ_NULL);
00247 
00248   /* */
00249   /* Get pointers to data */
00250   /* */
00251 
00252   /* Get the sub-vector we are reading from */
00253   assert(obj_data);
00254   state = (const struct RTOp_TOp_set_sub_vector_state_t*)obj_data;
00255 
00256   /* v (the vector to read from) */
00257   v_global_offset  = state->sub_vec.global_offset;
00258   v_sub_dim        = state->sub_vec.sub_dim;
00259   v_sub_nz         = state->sub_vec.sub_nz;
00260   v_val            = state->sub_vec.values;
00261   v_val_s          = state->sub_vec.values_stride;
00262   v_ind            = state->sub_vec.indices;
00263   v_ind_s          = state->sub_vec.indices_stride;
00264   v_l_off          = state->sub_vec.local_offset;
00265   v_sorted         = state->sub_vec.is_sorted;
00266 
00267   /* z (the vector to set) */
00268   z_global_offset  = targ_vecs[0].global_offset;
00269   z_sub_dim        = targ_vecs[0].sub_dim;
00270   z_val            = targ_vecs[0].values;
00271   z_val_s          = targ_vecs[0].values_stride;
00272 
00273   /* */
00274   /* Set part of the sub-vector for this chunk. */
00275   /* */
00276 
00277   if( v_global_offset + v_sub_dim < z_global_offset + 1
00278     || z_global_offset + z_sub_dim < v_global_offset + 1 )
00279     return 0; /* The sub-vector that we are setting does not overlap with this vector chunk! */
00280 
00281   if( v_sub_nz == 0 )
00282     return 0; /* The sparse sub-vector we are reading from is empty? */
00283 
00284   /* Get the number of elements that overlap */
00285   if( v_global_offset <= z_global_offset ) {
00286     if( v_global_offset + v_sub_dim >= z_global_offset + z_sub_dim )
00287       num_overlap = z_sub_dim;
00288     else
00289       num_overlap = (v_global_offset + v_sub_dim) - z_global_offset;
00290   }
00291   else {
00292     if( z_global_offset + z_sub_dim >= v_global_offset + v_sub_dim )
00293       num_overlap = v_sub_dim;
00294     else
00295       num_overlap = (z_global_offset + z_sub_dim) - v_global_offset;
00296   }
00297 
00298   /* Set the part of the sub-vector that overlaps */
00299   if( v_ind != NULL ) {
00300     /* Sparse elements */
00301     /* Set the overlapping elements to zero first. */
00302     if( v_global_offset >= z_global_offset )
00303       z_val += (v_global_offset - z_global_offset) * z_val_s;
00304     for( k = 0; k < num_overlap; ++k, z_val += z_val_s )
00305       *z_val = 0.0;
00306     /* Now set the sparse entries */
00307     z_val = targ_vecs[0].values;
00308     for( k = 0; k < v_sub_nz; ++k, v_val += v_val_s, v_ind += v_ind_s ) {
00309       i = v_global_offset + v_l_off + (*v_ind);
00310       if( z_global_offset < i && i <= z_global_offset + z_sub_dim )
00311         z_val[ z_val_s * (i - z_global_offset - 1) ] = *v_val;
00312     }
00313     /* ToDo: Implement a faster version for v sorted and eliminate the */
00314     /* if statement in the loop. */
00315   }
00316   else {
00317     /* Dense elemements */
00318     if( v_global_offset <= z_global_offset )
00319       v_val += (z_global_offset - v_global_offset) * v_val_s;
00320     else
00321       z_val += (v_global_offset - z_global_offset) * z_val_s;
00322     for( k = 0; k < num_overlap; ++k, v_val += v_val_s, z_val += z_val_s )
00323       *z_val = *v_val;
00324   }
00325 
00326   return 0; /* success? */
00327 }
00328 
00329 /* Public interface */
00330 
00331 const struct RTOp_RTOp_vtbl_t RTOp_TOp_set_sub_vector_vtbl =
00332 {
00333   &instance_obj_vtbl
00334   ,&RTOp_obj_null_vtbl /* use null type for reduction target object */
00335   ,"RTOp_TOp_set_sub_vector"
00336   ,NULL /* use default from reduct_vtbl */
00337   ,RTOp_TOp_set_sub_vector_apply_op
00338   ,NULL /* reduce_reduct_obj */
00339   ,NULL /* get_reduct_op */
00340 };
00341 
00342 int RTOp_TOp_set_sub_vector_construct(
00343   const struct RTOp_SparseSubVector* sub_vec, struct RTOp_RTOp* op )
00344 {
00345 #ifdef RTOp_DEBUG
00346   assert(sub_vec);
00347   assert(op);
00348 #endif
00349   op->vtbl     = &RTOp_TOp_set_sub_vector_vtbl;
00350   op->vtbl->obj_data_vtbl->obj_create(NULL,NULL,&op->obj_data);
00351   return RTOp_TOp_set_sub_vector_set_sub_vec(sub_vec,op);
00352 }
00353 
00354 int RTOp_TOp_set_sub_vector_set_sub_vec(
00355   const struct RTOp_SparseSubVector* sub_vec, struct RTOp_RTOp* op )
00356 {
00357   struct RTOp_TOp_set_sub_vector_state_t *state = NULL;
00358 #ifdef RTOp_DEBUG
00359   assert( op->obj_data );
00360 #endif
00361   state = (struct RTOp_TOp_set_sub_vector_state_t*)op->obj_data;
00362   if( state->owns_mem ) {
00363     free( (void*)state->sub_vec.values );
00364     free( (void*)state->sub_vec.indices );
00365   }
00366   state->sub_vec = *sub_vec;  /* We do not own the arrays values[] and indices[] */
00367   state->owns_mem = 0;
00368   return 0;
00369 }
00370 
00371 int RTOp_TOp_set_sub_vector_destroy( struct RTOp_RTOp* op )
00372 {
00373   op_free(NULL,NULL,&op->obj_data);
00374   op->vtbl = NULL;
00375   return 0;
00376 }

Generated on Thu Sep 18 12:33:39 2008 for RTOpPack: Extra C/C++ Code for Vector Reduction/Transformation Operators by doxygen 1.3.9.1