自学教程：C++ vset函数代码示例

/* * Ok, simulate a track-ball.  Project the points onto the virtual * trackball, then figure out the axis of rotation, which is the cross * product of P1 P2 and O P1 (O is the center of the ball, 0,0,0) * Note:  This is a deformed trackball-- is a trackball in the center, * but is deformed into a hyperbolic sheet of rotation away from the * center.  This particular function was chosen after trying out * several variations. * * It is assumed that the arguments to this routine are in the range * (-1.0 ... 1.0) */static voidtrackball ( float q[4], float p1x, float p1y, float p2x, float p2y ){    float a[3]; /* Axis of rotation */    float phi;  /* how much to rotate about axis */    float p1[3], p2[3], d[3];    float t;    if (p1x == p2x && p1y == p2y)    {        // Zero rotation        vzero(q);        q[3] = 1.0;        return;    }    // First, figure out z-coordinates for projection of P1 and P2 to    // deformed sphere    vset(p1,p1x,p1y,tb_project_to_sphere(TRACKBALLSIZE,p1x,p1y));    vset(p2,p2x,p2y,tb_project_to_sphere(TRACKBALLSIZE,p2x,p2y));    // Now, we want the cross product of P1 and P2    vcross(p2,p1,a);    // Figure out how much to rotate around that axis.    vsub(p1,p2,d);    t = vlength(d) / (2.0*TRACKBALLSIZE);    // Avoid problems with out-of-control values...    if (t > 1.0) t = 1.0;    if (t < -1.0) t = -1.0;    phi = 2.0 * asin(t);    axis_to_quat(a,phi,q);}

int DPMHC_K(struct str_DPMHC *ptr_DPMHC_data){    int i_K = ptr_DPMHC_data->i_K;    gsl_vector *v_u = ptr_DPMHC_data->v_u;    gsl_vector *v_v = ptr_DPMHC_data->v_v;    gsl_vector *v_w  = ptr_DPMHC_data->v_w;    gsl_matrix *m_DPtheta = ptr_DPMHC_data->m_DPtheta;    double d_DPalpha = ptr_DPMHC_data->d_DPalpha;    int K_tmp, K_new,j;    double a,v_j,w_j,csum,min_u;    //gsl_vector_view theta_j;  //int k_asset_number = P -> size1; /* number of assets in model */    K_tmp = i_K;    min_u = gsl_vector_min ( v_u );    a = 1.0 - min_u;    if( a == 1.0 )        printf("**********min_u = %g *************/n",min_u);    csum = 0.0;    j=0;    while ( csum <= a ){        /* check if new v_j,w_j and theta_j should be generated */        if( j >= K_tmp ){            v_j = gsl_ran_beta ( rng , 1.0, d_DPalpha );            vset( v_v, j, v_j);            w_j = v_j * (vget( v_w, j-1 )/vget(v_v,j-1))*(1.0-vget(v_v,j-1));            vset( v_w, j, w_j);        /* generate new mu, xi, tau from prior G_0 */            mset(m_DPtheta, j, 0,                ptr_DPMHC_data->d_m0 + gsl_ran_gaussian_ziggurat(rng, sqrt(ptr_DPMHC_data->d_s2m)));            mset(m_DPtheta, j, 1,                gsl_ran_gaussian_ziggurat(rng, ptr_DPMHC_data->d_A));            mset(m_DPtheta, j, 2,                gsl_ran_gamma(rng, 0.5, 0.5) );        }        csum += vget(v_w,j);        K_new = j + 1;        j++;    }    ptr_DPMHC_data->i_K = K_new;    return 0;}

object_t *vector_concat (object_t * a, object_t * b){  size_t al = VLENGTH (a), bl = VLENGTH (b);  object_t *c = c_vec (al + bl, NIL);  size_t i;  for (i = 0; i < al; i++)    vset (c, i, UPREF (vget (a, i)));  for (i = 0; i < bl; i++)    vset (c, i + al, UPREF (vget (b, i)));  return c;}

void init()                                         // --- Initialization routine{    double tension = c_Mpay;    double fcoef = 0;    double mt = 0;    curtime = 0;                                    // Set current time to 0    numpoints = NUMPOINTS;                          // Set variable to macro    points = smalloc((numpoints + 1) * sizeof(point_t)); // Allocate points (FIXME)    for(int i = 0; i < numpoints; i++) {            // Loop over points        vset(&points[i].pos, c_Clen * i / NUMPOINTS + c_Re, 0, 0);        vset(&points[i].v, 0, 0, 0);        points[i].crosssect = c_Fs * tension / c_Slim;        points[i].k = c_Ecnt / (c_Clen / NUMPOINTS) * points[i].crosssect;        points[i].initlen = (c_Clen / NUMPOINTS - tension / points[i].k);        points[i].frict = c_Fr_cnt;        mt += points[i].m = points[i].initlen * points[i].crosssect * c_Dcnt;        fcoef = c_Ge / (points[i].pos.x * points[i].pos.x) - points[i].pos.x * c_Vang * c_Vang;        tension += fcoef * points[i].m;    }    points[numpoints - 1].m = -tension / fcoef * 1; // FIXME    if(0) {        for(int i = 0; i < numpoints; i++) {           printf("Item number: %d/n", i + 1);           printf("Position: %f, %f, %f/n", points[i].pos.x, points[i].pos.y, points[i].pos.z);           printf("Velocity: %f, %f, %f/n", points[i].v.x, points[i].v.y, points[i].v.z);           printf("Force: %f, %f, %f/n", points[i].f.x, points[i].f.y, points[i].f.z);           printf("Mass: %f/n", points[i].m);           printf("Spring Constant: %f/n", points[i].k);           printf("Friction Coefficient: %f/n", points[i].frict);           printf("Initial Length: %f/n", points[i].initlen);           printf("Cross Sectional Area: %f/n", points[i].crosssect);        }    } else if(1) {        for(int i = 0; i < numpoints; i++) {           //printf("[%d, ", i + 1);           printf("[(%e, %e, %e), ", points[i].pos.x, points[i].pos.y, points[i].pos.z);           printf("(%e, %e, %e), ", points[i].v.x, points[i].v.y, points[i].v.z);           printf("(%e, %e, %e), ", points[i].f.x, points[i].f.y, points[i].f.z);           printf("%e, ", points[i].m);           printf("%e, ", points[i].k);           printf("%e, ", points[i].frict);           printf("%e, ", points[i].initlen);           printf("%e]/n", points[i].crosssect);        }    }}

bool SlideNavmesh::mouseDown(const float x, const float y){	if (!m_expanded)	{		if (hitCorner(x,y))		{			startDrag(1, x,y, m_pos);			return true;		}	}	else	{		int bidx = hitButtons(x-m_pos[0],y-m_pos[1]);		if (bidx != -1)		{			return true;		}		else if (hitCorner(x,y))		{			startDrag(1, x,y, m_pos);			return true;		}		else if (hitArea(x,y))		{			const float lx = x - (m_pos[0]+PADDING_SIZE);			const float ly = y - (m_pos[1]+PADDING_SIZE);			float pos[2] = {lx,ly};			float nearest[2] = {lx,ly};			if (m_scene.nav)				navmeshFindNearestTri(m_scene.nav, pos, nearest);						if (SDL_GetModState() & KMOD_SHIFT)			{				agentMoveAndAdjustCorridor(&m_scene.agents[0], nearest, m_scene.nav);				vcpy(m_scene.agents[0].oldpos, m_scene.agents[0].pos);				vset(m_scene.agents[0].corner, FLT_MAX,FLT_MAX);			}			else			{				vcpy(m_scene.agents[0].target, nearest);				vcpy(m_scene.agents[0].oldpos, m_scene.agents[0].pos);				agentFindPath(&m_scene.agents[0], m_scene.nav);				vset(m_scene.agents[0].corner, FLT_MAX,FLT_MAX);			}			return true;		}	}	return false;}

int pod_experiment(char* observed_data, char* heldout_data,                   char* model_root, char* out){    corpus *obs, *heldout;    llna_model *model;    llna_var_param *var;    int i;    gsl_vector *log_lhood, *e_theta;    doc obs_doc, heldout_doc;    char string[100];    double total_lhood = 0, total_words = 0, l;    FILE* e_theta_file = fopen("/Users/blei/llna050_e_theta.txt", "w");    // load model and data    obs = read_data(observed_data);    heldout = read_data(heldout_data);    assert(obs->ndocs == heldout->ndocs);    model = read_llna_model(model_root);    // run experiment    init_temp_vectors(model->k-1); // !!! hacky    log_lhood = gsl_vector_alloc(obs->ndocs + 1);    e_theta = gsl_vector_alloc(model->k);    for (i = 0; i < obs->ndocs; i++)    {        // get observed and heldout documents        obs_doc = obs->docs[i];        heldout_doc = heldout->docs[i];        // compute variational distribution        var = new_llna_var_param(obs_doc.nterms, model->k);        init_var_unif(var, &obs_doc, model);        var_inference(var, &obs_doc, model);        expected_theta(var, &obs_doc, model, e_theta);        vfprint(e_theta, e_theta_file);        // approximate inference of held out data        l = log_mult_prob(&heldout_doc, e_theta, model->log_beta);        vset(log_lhood, i, l);        total_words += heldout_doc.total;        total_lhood += l;        printf("hid doc %d    log_lhood %5.5f/n", i, vget(log_lhood, i));        // save results?        free_llna_var_param(var);    }    vset(log_lhood, obs->ndocs, exp(-total_lhood/total_words));    printf("perplexity : %5.10f", exp(-total_lhood/total_words));    sprintf(string, "%s-pod-llna.dat", out);    printf_vector(string, log_lhood);    return(0);}

int DPMHC_v_smplr(struct str_DPMHC *ptr_DPMHC_data){    gsl_vector *v_v = ptr_DPMHC_data->v_v;    gsl_vector *v_w = ptr_DPMHC_data->v_w;    gsl_vector_int *vi_S = ptr_DPMHC_data->vi_S;    int i_K = ptr_DPMHC_data->i_K;    double d_DPalpha = ptr_DPMHC_data->d_DPalpha;    size_t i_T = vi_S->size;    int i,j,i_si;    double d_vj,d_prod;    //   printf("inside sample_v K=%d/n",K);    gsl_vector *v_a = gsl_vector_alloc ( (size_t) i_K );    gsl_vector *v_b = gsl_vector_alloc ( (size_t) i_K );    gsl_vector_set_all ( v_a, 1.0 );    gsl_vector_set_all ( v_b, d_DPalpha);    for(i=0;i<i_T;i++){        i_si = vget_int(vi_S,i);        (v_a->data[ i_si ]) += 1.0;        for(j = 0; j < i_si;j++){            (v_b->data[ j ]) += 1.0;        }    }   /* pvec(a); */   /* pvec(b); */    /* take draws and form w */    for(j=0;j<i_K;j++){        d_vj = gsl_ran_beta ( rng , vget(v_a,j) , vget(v_b,j) );        vset( v_v, j, d_vj );        /* w_j */        if( j == 0 ){            vset(v_w,j, d_vj );            d_prod = (1.0-d_vj);        }else{            vset(v_w,j, d_vj*d_prod );            d_prod *= (1.0-d_vj);        }    }    gsl_vector_free (v_a);    gsl_vector_free (v_b);    return 0;}

struct summary * summarise_vec( VEC v){	assert(NULL!=v);	VEC quant = create_vec(5);	vset(quant,0,0.); vset(quant,1,0.25); vset(quant,2,0.5); vset(quant,3,0.75); vset(quant,4,1.);	struct summary * s = malloc(sizeof(struct summary));	s->mean = mean(v);	s->var = variance(v);	s->quantiles = quantiles(v,quant);	s->mad = mad(v);	s->data = v;	return s;}

double compute_lda_lhood(lda_post* p) {  int k, n;  int K = p->model->ntopics, N = p->doc->nterms;  double gamma_sum = sum(p->gamma);  double lhood =    gsl_sf_lngamma(sum(p->model->alpha)) -    gsl_sf_lngamma(gamma_sum);  vset(p->lhood, K, lhood);    double influence_term = 0.0;  double digsum = gsl_sf_psi(gamma_sum);  for (k = 0; k < K; k++) {    if (p->doc_weight != NULL) {      //	  outlog("doc weight size: %d", p->doc_weight->size);      assert (K == p->doc_weight->size);      double influence_topic = gsl_vector_get(p->doc_weight, k);      if (FLAGS_model == "dim"	  || FLAGS_model == "fixed") {	influence_term = - ((influence_topic * influence_topic			     + FLAGS_sigma_l * FLAGS_sigma_l)			    / 2.0 / (FLAGS_sigma_d * FLAGS_sigma_d));	// Note that these cancel with the entropy.	//     - (log(2 * PI) + log(FLAGS_sigma_d)) / 2.0);      }    }    double e_log_theta_k = gsl_sf_psi(vget(p->gamma, k)) - digsum;    double lhood_term =      (vget(p->model->alpha, k)-vget(p->gamma, k)) * e_log_theta_k +      gsl_sf_lngamma(vget(p->gamma, k)) -      gsl_sf_lngamma(vget(p->model->alpha, k));        for (n = 0; n < N; n++) {      if (mget(p->phi, n, k) > 0) {	lhood_term +=	  p->doc->count[n]*	  mget(p->phi, n, k) *	  (e_log_theta_k	   + mget(p->model->topics, p->doc->word[n], k)	   - mget(p->log_phi, n, k));      }    }    vset(p->lhood, k, lhood_term);    lhood += lhood_term;    lhood += influence_term;  }    return(lhood);}

void simul(){    for(int i = 0; i < numpoints; i++) { vset(&points[i].f, 0, 0, 0); }    for(int i = 0; i < numpoints - 1; i++) {        vector_t vdiff;        GLdouble axialv, frict, elast;        vector_t diff = vsub(&points[i+1].pos, &points[i].pos);        GLdouble l = vnorm(&diff);        vector_t f;        elast = points[i].k * (1 - points[i].initlen / l);        if(elast < 0) { elast = 0; }        points[i].tension = elast * l;        vdiff = vsub(&points[i+1].v, &points[i].v);        axialv = vdot(&vdiff, &diff) / l;        frict = axialv / l * points[i].frict;        f = vscale(elast + frict, &diff);        vaddto(&points[i].f, &f);        vsubfrom(&points[i+1].f, &f);    }    for(int i = 0; i < numpoints; i++) {        double loading = points[i].tension / points[i].crosssect / c_Slim;        if(loading > 1 && curtime > 0) {            points[i].k = 0;            points[i].frict = 0;        }    }    for(int i = 0; i < numpoints; i++) { pointdynamics(&points[i]); }    {        points[0].pos.x = c_Re * cos(0 / 180.0 * M_PI);        points[0].pos.y = 0;        points[0].pos.z = c_Re * sin(0 / 180.0 * M_PI);        vset(&points[0].v, 0, 0, 0);        if(curtime > 1) {            points[BREAKPOINT].k = 0;            points[BREAKPOINT].frict = 0;        }    }    curtime += TIMESTEP;}

int DPMHC_u_smplr(struct str_DPMHC *ptr_DPMHC_data){    gsl_vector_int *vi_S = ptr_DPMHC_data->vi_S;    gsl_vector *v_w = ptr_DPMHC_data->v_w;    gsl_vector *v_u = ptr_DPMHC_data->v_u;    size_t i_T = vi_S->size;    int i,i_si;    double d_ui,d_w_si;    for(i=0;i<i_T;i++){        i_si = vget_int( vi_S, i);        d_w_si = vget( v_w, i_si);        dw:            d_ui = d_w_si * gsl_rng_uniform (rng);            if( d_ui == 0.0 ){            printf("/n/n ******** u_i = 0.0 /n");            goto dw;            }        vset( v_u, i, d_ui);    }    return 0;}

float *calc_arc(float center[3],double radius,int divisions, int start, int end){	float *points = (float *)mem_malloc(sizeof(float)*(divisions*3));	int i;	int j=0;	int angle = end - start;	float phi = deg_to_rad(angle);	double delta = (double)((double)(phi ) / (divisions -1 ));	double a = 0;	for (i=0;i<divisions;i++)	{		float r[3];		float result[3];		vset(r,cos(a),sin(a),0);		a += delta;		vmul(r,(float)radius);		vadd(result,center,r);		points[j]=result[0];		points[j+1]=result[1];		points[j+2]=result[2];		j+=3;	}	return points;}

void vcross( triple* x, const triple* y ){  triple t;  vset( &t, x );  x->x = ( t.y * y->z - t.z * y->y );  x->y = ( t.z * y->x - t.x * y->z );  x->z = ( t.x * y->y - t.y * y->x );}

inline void pointdynamics(point_t *curpt)           // Calculate point movement{    GLdouble r = vnorm(&curpt->pos);    vector_t grav;    vector_t rotf;    vector_t coriolis;    if(r > c_Re) {                                  // Points outside Earth        grav = vscale(-c_Ge * curpt->m / pow(r,3), &curpt->pos);        vaddto(&curpt->f, &grav);    } else {                                        // Take care of points that        vector_t scaled;                            // are inside Earth        curpt->pos = vscale(c_Re / vnorm(&curpt->pos), &curpt->pos);        scaled = vscale(-vdot(&curpt->pos, &curpt->v) / c_Re / c_Re, &curpt->pos);        vaddto(&curpt->v, &scaled);        curpt->v = vscale(0, &curpt->v);    }    rotf = vscale(c_Vang * c_Vang * curpt->m, &curpt->pos);    rotf.z = 0;    vset(&coriolis, 2 * curpt->v.y * curpt->m * c_Vang, -2 * curpt->v.x * curpt->m * c_Vang, 0);    vaddto(&curpt->f, &rotf);    vaddto(&curpt->f, &coriolis);    vaddscaled(&curpt->v, TIMESTEP / curpt->m, &curpt->f);    vaddscaled(&curpt->pos, TIMESTEP, &curpt->v);    if(curtime <= 0) { curpt->v = vscale(.99, &curpt->v); }}

double lda_m_step(lda* model, lda_suff_stats* ss) {    int k, w;    double lhood = 0;    for (k = 0; k < model->ntopics; k++)    {        gsl_vector ss_k = gsl_matrix_column(ss->topics_ss, k).vector;        gsl_vector log_p = gsl_matrix_column(model->topics, k).vector;        if (LDA_USE_VAR_BAYES == 0)        {            gsl_blas_dcopy(&ss_k, &log_p);            normalize(&log_p);            vct_log(&log_p);        }        else        {            double digsum = sum(&ss_k)+model->nterms*LDA_TOPIC_DIR_PARAM;            digsum = gsl_sf_psi(digsum);            double param_sum = 0;            for (w = 0; w < model->nterms; w++)            {                double param = vget(&ss_k, w) + LDA_TOPIC_DIR_PARAM;                param_sum += param;                double elogprob = gsl_sf_psi(param) - digsum;                vset(&log_p, w, elogprob);                lhood += (LDA_TOPIC_DIR_PARAM - param) * elogprob + gsl_sf_lngamma(param);            }            lhood -= gsl_sf_lngamma(param_sum);        }    }    return(lhood);}

void update_phi(int doc_number, int time,		lda_post* p, lda_seq* var,		gsl_matrix* g) {    int i, k, n, K = p->model->ntopics, N = p->doc->nterms;    double dig[p->model->ntopics];    for (k = 0; k < K; k++) {      dig[k] = gsl_sf_psi(vget(p->gamma, k));    }    for (n = 0; n < N; n++) {      // compute log phi up to a constant      int w = p->doc->word[n];      for (k = 0; k < K; k++) {	mset(p->log_phi, n, k,	     dig[k] + mget(p->model->topics, w, k));      }      // normalize in log space      gsl_vector log_phi_row = gsl_matrix_row(p->log_phi, n).vector;      gsl_vector phi_row = gsl_matrix_row(p->phi, n).vector;      log_normalize(&log_phi_row);      for (i = 0; i < K; i++) {	vset(&phi_row, i, exp(vget(&log_phi_row, i)));      }    }}

void SlideNavmesh::update(const float dt){	if (!m_update && !m_step)		return;	m_step = false;	const float maxSpeed = 1.0f;	NavmeshAgent* agent = &m_scene.agents[0];		// Find next corner to steer to.	// Smooth corner finding does a little bit of magic to calculate spline	// like curve (or first tangent) based on current position and movement direction	// next couple of corners.	float corner[2],dir[2];	int last = 1;	vsub(dir, agent->pos, agent->oldpos); // This delta handles wall-hugging better than using current velocity.	vnorm(dir);	vcpy(corner, agent->pos);	if (m_moveMode == AGENTMOVE_SMOOTH || m_moveMode == AGENTMOVE_DRUNK)		last = agentFindNextCornerSmooth(agent, dir, m_scene.nav, corner);	else		last = agentFindNextCorner(agent, m_scene.nav, corner);			if (last && vdist(agent->pos, corner) < 0.02f)	{		// Reached goal		vcpy(agent->oldpos, agent->pos);		vset(agent->dvel, 0,0);		vcpy(agent->vel, agent->dvel);		return;	}	vsub(agent->dvel, corner, agent->pos);	// Apply style	if (m_moveMode == AGENTMOVE_DRUNK)	{		agent->t += dt*4;		float amp = cosf(agent->t)*0.25f;		float nx = -agent->dvel[1];		float ny = agent->dvel[0];		agent->dvel[0] += nx * amp;		agent->dvel[1] += ny * amp;	}		// Limit desired velocity to max speed.	const float distToTarget = vdist(agent->pos,agent->target);	const float clampedSpeed = maxSpeed * min(1.0f, distToTarget/agent->rad);	vsetlen(agent->dvel, clampedSpeed);	vcpy(agent->vel, agent->dvel);	// Move agent	vscale(agent->delta, agent->vel, dt);	float npos[2];	vadd(npos, agent->pos, agent->delta);	agentMoveAndAdjustCorridor(&m_scene.agents[0], npos, m_scene.nav);}

void center(gsl_vector* v){    int size = v->size;    double mean = sum(v)/size;    int i;    for (i = 0; i < size; i++)        vset(v, i, vget(v,i)-mean);}

void normalize(gsl_vector* v){    int size = v->size;    double sum_v = sum(v);    int i;    for (i = 0; i < size; i++)        vset(v, i, vget(v,i)/sum_v);}

/* * Ok, simulate a track-ball.  Project the points onto the virtual * trackball, then figure out the axis of rotation, which is the cross * product of P1 P2 and O P1 (O is the center of the ball, 0,0,0) * Note:  This is a deformed trackball-- is a trackball in the center, * but is deformed into a hyperbolic sheet of rotation away from the * center.  This particular function was chosen after trying out * several variations. * * It is assumed that the arguments to this routine are in the range * (-1.0 ... 1.0) */void trackball( double q[4], double p1x, double p1y, double p2x, double p2y ){    double a[3]; /* Axis of rotation */    double phi;  /* how much to rotate about axis */    double p1[3], p2[3], d[3];    double t;    if( p1x == p2x && p1y == p2y )    {        /* Zero rotation */        vzero( q );        q[3] = 1.0;        return;    }    /*     * First, figure out z-coordinates for projection of P1 and P2 to     * deformed sphere     */    vset( p1, p1x, p1y, tb_project_to_sphere( TRACKBALLSIZE, p1x, p1y ) );    vset( p2, p2x, p2y, tb_project_to_sphere( TRACKBALLSIZE, p2x, p2y ) );    /*     *  Now, we want the cross product of P1 and P2     */    vcross(p2,p1,a);    /*     *  Figure out how much to rotate around that axis.     */    vsub( p1, p2, d );    t = vlength( d ) / (2.0f * TRACKBALLSIZE);    /*     * Avoid problems with out-of-control values...     */    if( t > 1.0 )        t = 1.0;    if( t < -1.0 )        t = -1.0;    phi = 2.0f * (double) asin( t );    axis_to_quat( a, phi, q );}

object_t *vset_check (object_t * vo, object_t * io, object_t * val){  int i = into2int (io);  vector_t *v = OVAL (vo);  if (i < 0 || i >= (int) v->len)    THROW (out_of_bounds, UPREF (io));  vset (vo, i, UPREF (val));  return UPREF (val);}

/* * Scan attribute sets 'pa' and 'pb' for common vendors and * call their file_attributes_match() methods to find incompatibilities. * * If the vendor sets in 'pa' and 'pb' do not match or if only one set * exists (the other being NULL) this is treated as a mismatch. * * RETURNS: zero if no incompatibility is found, a value < 0 otherwise. * * NOTE:    it is legal to pass 'pa==NULL', 'pb==NULL' in which case *          the routine returns 0 (SUCCESS). */intpmelf_match_attribute_set(Pmelf_attribute_set *pa, Pmelf_attribute_set *pb){int                    i,j;Pmelf_attribute_vendor *pv;int                    vendor_set_a, vendor_set_b;int                    rval = 0;	if ( !pa && !pb )		return 0;	if ( ! (pa && pb) )		return -1;	if ( (vendor_set_a = vset(pa)) < 0 ) {		return -1;	}	if ( (vendor_set_b = vset(pb)) < 0 ) {		return -1;	}	if ( vendor_set_a != vendor_set_b ) {		PMELF_PRINTF(pmelf_err, PMELF_PRE"pmelf_match_attribute_set(): vendor set mismatch (objects %s, %s)/n", pa->obj_name, pb->obj_name);		return -2;	}	for ( i=0; i<ATTR_MAX_VENDORS; i++ ) {		if ( pa->attributes[i].file_attributes ) {			pv = pa->attributes[i].file_attributes->pv;			for ( j=0; j<ATTR_MAX_VENDORS; j++ ) {				if ( pb->attributes[j].file_attributes ) {					if ( pv == pb->attributes[j].file_attributes->pv ) {						if ( ! pv->file_attributes_match ) {							PMELF_PRINTF(pmelf_err, PMELF_PRE"pmelf_match_attribute_set(): vendor %s has no 'match' handler/n", pv->name);							return -1;						}						rval |= pv->file_attributes_match( pa->attributes[i].file_attributes, pb->attributes[j].file_attributes );					}				}			}		}	}	return rval;}

/* * normalize a vector in log space * * x_i = log(a_i) * v = log(a_1 + ... + a_k) * x_i = x_i - v * */void log_normalize(gsl_vector* x) {	double v = vget(x, 0);	for (unsigned int i = 1; i < x->size; i++) {		v = log_sum(v, vget(x, i));	}		for (unsigned int i = 0; i < x->size; i++) {		vset(x, i, vget(x,i)-v);	}}

VEC branchlengths_from_tree ( const TREE * tree){	assert(NULL!=tree);		VEC branlength = create_vec(tree->n_br);	assert(NULL!=branlength);	for (unsigned int bran=0 ; bran<tree->n_br ; bran++){		vset(branlength,bran,tree->branches[bran]->blength[0]);	}	return branlength;}

void normalize(gsl_vector* x) {	double v = 0;		for (unsigned int i = 0; i < x->size; i++) {		v += vget(x, i);	}		for (unsigned int i = 0; i < x->size; i++) {		vset(x, i, vget(x, i) / v);	}}

object_t *vector_sub (object_t * vo, int start, int end){  vector_t *v = OVAL (vo);  if (end == -1)    end = v->len - 1;  object_t *newv = c_vec (1 + end - start, NIL);  int i;  for (i = start; i <= end; i++)    vset (newv, i - start, UPREF (vget (vo, i)));  return newv;}

void init_lda_post(lda_post* p) {    int k, n, K = p->model->ntopics, N = p->doc->nterms;    for (k = 0; k < K; k++)    {        vset(p->gamma, k,             vget(p->model->alpha,k) + ((double) p->doc->total)/K);        for (n = 0; n < N; n++)            mset(p->phi, n, k, 1.0/K);    }    p->doc_weight = NULL;}

void initialize_lda_ss_from_random(corpus_t* data, lda_suff_stats* ss) {    int k, n;    gsl_rng * r = new_random_number_generator();    for (k = 0; k < ss->topics_ss->size2; k++)    {        gsl_vector topic = gsl_matrix_column(ss->topics_ss, k).vector;        gsl_vector_set_all(&topic, 0);        for (n = 0; n < topic.size; n++)        {	  vset(&topic, n, gsl_rng_uniform(r) + 0.5 / data->ndocs + 4.0);        }    }}