C++ IVector类代码示例

void zad1() {
    IVector *v1 = (Vector::parseSimple("2 3 -4"))->add(Vector::parseSimple("-1 4 -3"));
    printf("v1 = %s\n", v1->toString().c_str());
    
    IVector *s = v1->nAdd(Vector::parseSimple("-1 4 -3"));
    printf("s = %s\n", s->toString().c_str());
    
    IVector *v2 = v1->nVectorProduct(Vector::parseSimple("2 2 4"));
    printf("v2 = %s\n", v2->toString().c_str());
    
    IVector *v3 = v2->normalize();
    printf("v3 = %s\n", v3->toString().c_str());
    
    IVector *v4 = v2->scalarMultiply(-1);
    printf("v4 = %s\n", v4->toString().c_str());
    
    IMatrix *m1 = (Matrix::parseSimple("1 2 3 | 2 1 3 | 4 5 1"))
    ->add(Matrix::parseSimple("-1 2 -3 | 5 -2 7 | -4 -1 3"));
    printf("m1 = %s\n", m1->toString().c_str());
    
    IMatrix *m2 = (Matrix::parseSimple("1 2 3 | 2 1 3 | 4 5 1"))
    ->nMultiply(Matrix::parseSimple("-1 2 -3 | 5 -2 7 | -4 -1 3")->nTranspose(false));
    printf("m2 = %s\n", m2->toString().c_str());
    
    IMatrix *m3 = (Matrix::parseSimple("-24 18 5 | 20 -15 -4 | -5 4 1")->nInvert())
    ->nMultiply(Matrix::parseSimple("1 2 3 | 0 1 4 | 5 6 0")->nInvert());
    printf("m3 = %s\n", m3->toString().c_str());
}

int ReferenceFrame2D::GetNode(Vector2D p_loc) const
{
	double tol = 1e-8*GetMBS()->CharacteristicLength();
/*
	for (int i = 1; i <=nodes.Length(); i++)
	{
		if (Dist(p_loc,nodes(i)->Pos2D()) <= tol) return i;
	}
	GetMBS()->UO() << "ERROR: node not found:" << p_loc << " !!!\n";
	return 0;
*/

	IVector items;
	Vector3D p(p_loc.X(),p_loc.Y(),0.);

	Box3D box(p, p);
	box.Increase(tol);
	searchtree.GetItemsInBox(box,items);
	for (int i = 1; i <= items.Length(); i++)
	{
		Vector3D p2(nodes(items(i))->Pos().X(), nodes(items(i))->Pos().Y(), 0.);
		if (Dist(p,p2) <= tol) return items(i);
	}
	GetMBS()->UO() << "ERROR: node not found:" << p_loc << " !!!\n";
	return 0;
}

//The function recursively builds up the transformation of the pThing 
//until it's in room-space coordinates
HRESULT C2DThingCoordTransformer::BuildTransformation(IThing * pThing, IThing * pParentThing)
{
	IVector * pVector = NULL;
	HRESULT hr = S_OK;
	IThing* pTmpParentThing = NULL;
	D3DRMMATRIX4D tmpMatrix;
	float flScaleX, flScaleY, flScaleZ;
	float flThingX, flThingY, flThingZ, // Position of current Thing having cells edited
	      flThingDirX, flThingDirY, flThingDirZ; // Orientation of current Thing having cells edited

	if (!pThing || !pParentThing)
		goto EXIT_FAIL;

	// Store the position Vector of the Thing
	hr = pThing->get_ObjectProperty(bstrPosition, (IObjectProperty **) &pVector);
	if( FAILED(hr) || !pVector) goto EXIT_FAIL;
	hr = pVector->get(&flThingX, &flThingY, &flThingZ);
	if( FAILED(hr) ) goto EXIT_FAIL;
	SAFERELEASE(pVector);

	// Get the Scale Vector of the Thing
	hr = pThing->get_ObjectProperty(bstrScale, (IObjectProperty **) &pVector);
	if( FAILED(hr) || !pVector) goto EXIT_FAIL;
	hr = pVector->get(&flScaleX, &flScaleY, &flScaleZ);
	if( FAILED(hr) ) goto EXIT_FAIL;
	SAFERELEASE(pVector);

	// Store the Direction Vector of the Thing
	hr = pThing->get_ObjectProperty(bstrDirection, (IObjectProperty **) &pVector);
	if( FAILED(hr) ) goto EXIT_FAIL;
	hr = pVector->get(&flThingDirX, &flThingDirY, &flThingDirZ);
	if( FAILED(hr) ) goto EXIT_FAIL;
	SAFERELEASE(pVector);

	IdentityMatrix(&tmpMatrix);
	ScaleMatrix(&tmpMatrix, flScaleX, flScaleY, flScaleZ);
	RotateMatrix(&tmpMatrix, flThingDirX, flThingDirY, flThingDirZ);
	TranslateMatrix(&tmpMatrix, flThingX, flThingY, flThingZ);
	PostMultiplyMatrix(&m_d3dMatrix, &tmpMatrix);

	IdentityMatrix(&tmpMatrix);
	TranslateMatrix(&tmpMatrix, -flThingX, -flThingY, -flThingZ);
	InverseRotateMatrix(&tmpMatrix, flThingDirX, flThingDirY, flThingDirZ);
	ScaleMatrix(&tmpMatrix, 1.0f / flScaleX, 1.0f / flScaleY, 1.0f / flScaleZ);
	PostMultiplyMatrix(&m_d3dInverseMatrix, &tmpMatrix);

	hr = pParentThing->get_Container(&pTmpParentThing);
	if ( SUCCEEDED(hr) && pTmpParentThing)
	{
		BuildTransformation(pParentThing, pTmpParentThing);
	}


EXIT_FAIL:
	SAFERELEASE(pTmpParentThing);
	SAFERELEASE(pVector);
	return hr;
}

void verifyExistenceOfPeriodicOrbit(IPoincareMap& pm, IVector X, int period)
{
  IVector center = midVector(X);
  interval returnTime;

  // Center is 2-dimensional. Embed it into Poincare section, i.e. add first coordinate x=0.
  // Define a tripleton representation of the center of X.
  C0HOTripletonSet s1({interval(0.),center[0],center[1]});

  // Compute iteration of Poincare map at the center (3-dim object is returned).
  IVector y = pm(s1,returnTime,period);

  // Project it onto 2-dim Poincare section, first coordinate is ignored.
  IVector imCenter(2,y.begin()+1);


  // Derivative of PM on the set X.
  // Define doubleton representation of the first order variational equations.
  C1HORect2Set s2({interval(0.),X[0],X[1]});

  // The matrix monodromyMatrix will store derivatives of the FLOW not Poincare map.
  IMatrix monodromyMatrix(3,3);
  y = pm(s2,monodromyMatrix,returnTime,period);

  // This member function recomputes derivatives of the flow into derivatives of Poincare map
  IMatrix DP = pm.computeDP(y,monodromyMatrix);

  // We extract a 2x2 slice from 3x3 DP matrix and subtract identity.
  IMatrix DP_Minus_Id(2,2);
  DP_Minus_Id[0][0] = DP[1][1] - 1.;
  DP_Minus_Id[0][1] = DP[1][2];
  DP_Minus_Id[1][0] = DP[2][1];
  DP_Minus_Id[1][1] = DP[2][2] - 1.;

  // Compute interval Newton operator.
  IVector N = center - capd::matrixAlgorithms::gauss(DP_Minus_Id,imCenter-center);

  // Verification if N is a subset of X
  cout << "\n---------------------------------------------------\n\nN = " << N << endl;
  cout << "X = " << X << endl;
  cout << "Return time: " << returnTime << endl;
  if(subsetInterior(N,X))
    cout << "the existence of period " << period << " orbit verified\n";
  else
  {
    cout << "N is not a subset of X\n\n";
    cout << "diam(N)=" << diam(N) << endl;
    cout << "diam(X)=" << diam(X) << endl;
    cout << "N-X" << N-X << endl;
  }
}

 void izracunajKoeficijente() {
     for(int i = 0; i < (int)faces.size(); ++i) {
         int i0 = faces[i].indexes[0];
         int i1 = faces[i].indexes[1];
         int i2 = faces[i].indexes[2];
         IVector *v1 = new Vector(new double[3]{vertices[i1].x - vertices[i0].x, vertices[i1].y - vertices[i0].y, vertices[i1].z - vertices[i0].z }, 3);
         IVector *v2 = new Vector(new double[3]{vertices[i2].x - vertices[i0].x, vertices[i2].y - vertices[i0].y, vertices[i2].z - vertices[i0].z }, 3);
         IVector *n = v1->nVectorProduct(v2);
         faces[i].a = n->get(0);
         faces[i].b = n->get(1);
         faces[i].c = n->get(2);
         faces[i].d = -faces[i].a * vertices[i1].x - faces[i].b * vertices[i1].y - faces[i].c * vertices[i1].z;
     }
 }

    static auto prepare_it(IIterator ifirst, IIterator ilast, EIterator efirst, EIterator elast, IVector& ivec, EVector& evec) {
        std::copy(ifirst, ilast, std::back_inserter(ivec));

        auto input_first = ivec.begin();
        auto input_last  = ivec.end();

        if (Denoising) {
            std::copy(efirst, elast, std::back_inserter(evec));

            auto expected_first = evec.begin();
            auto expected_last = evec.end();
            return std::make_tuple(input_first, input_last, expected_first, expected_last);
        } else {
            return std::make_tuple(input_first, input_last, input_first, input_last);
        }
    }

ode_solver::ODE_result ode_solver::simple_ODE_backward(IVector & X_0, IVector const & X_t, interval const & T,
                                                       IVector const & inv, vector<IFunction> & funcs) {
    // X_0 = X_0 \cup (X_t - + (d/dt Inv) * T)
    for (int i = 0; i < X_0.dimension(); i++) {
        interval & x_0 = X_0[i];
        interval const & x_t = X_t[i];
        IFunction & dxdt = funcs[i];
        for (Enode * par : m_pars) {
            double lb = get_lb(par);
            double ub = get_ub(par);
            string name = par->getCar()->getName();
            dxdt.setParameter(name, interval(lb, ub));
        }
        try {
            interval const new_x_0 = x_t - dxdt(inv) * T;
            if (!intersection(new_x_0, x_0, x_0)) {
                DREAL_LOG_INFO << "Simple_ODE: no intersection for X_0";
                return ODE_result::UNSAT;
            }
        } catch (exception& e) {
            DREAL_LOG_INFO << "Exception in Simple_ODE: " << e.what();
        }
    }
    // update
    IVector_to_varlist(X_0, m_0_vars);
    return ODE_result::SAT;
}

void ode_solver::IVector_to_varlist(IVector const & v, vector<Enode*> & vars) {
    for (auto i = 0; i < v.dimension(); i++) {
        double lb = get_lb(vars[i]);
        double ub = get_ub(vars[i]);
        if (lb < v[i].leftBound())
            set_lb(vars[i], v[i].leftBound());
        if (ub > v[i].rightBound())
            set_ub(vars[i], v[i].rightBound());
    }
}

void printSample(const double &t, const IVector &x)
{
	cout.setf(ios::scientific);
  cout.setf(ios::showpos);
  cout.precision(10);
  cout << t << " ";
  for(int i=0; i<x.dimension(); ++i)
    cout << x[i].leftBound() << " " << x[i].rightBound() << " ";
  cout << endl; 
}

// Take an intersection of v and inv.
// If there is no intersection, return false.
bool ode_solver::check_invariant(IVector & v, IVector const & inv) {
    if (!intersection(v, inv, v)) {
        DREAL_LOG_INFO << "invariant violated!";
        for (auto i = 0; i < v.dimension(); i++) {
            if (v[i].leftBound() < inv[i].leftBound() || v[i].rightBound() > inv[i].rightBound()) {
                DREAL_LOG_INFO << "inv[" << i << "] = " << inv[i];
                DREAL_LOG_INFO << "  v[" << i << "] = " <<   v[i];
            }
        }
        return false;
    }
    return true;
}

int ReferenceFrame2D::AddNode(Node* n)
{
	n->SetAuxNode();

	double tol = 1e-8*GetMBS()->CharacteristicLength();

	//standard:
	/*
	for (int i = 1; i <=nodes.Length(); i++)
	{
		if (Dist(n->Pos(),nodes(i)->Pos()) <= tol) return i;
	}
	Node* nc = n->GetCopy();
	return nodes.Add(nc);
	*/

	//optimized:

	IVector items;

	Box3D box(n->Pos(),n->Pos());
	box.Increase(tol); //disabled in example for witteven adams comparison Milano conference 2007

	searchtree.GetItemsInBox(box,items);
	for (int i = 1; i <= items.Length(); i++)
	{
		//if (Dist(n->Pos(),nodes(items(i))->Pos()) <= tol) return items(i);
		if (Dist(n->Pos(),nodes(items(i))->Pos()) <= tol && (n->GetBodyInd() == nodes(items(i))->GetBodyInd()) ) return items(i); //$ AD FENodes have domain
	}
	Node* nc = n->GetCopy();
	int index = nodes.Add(nc);
	nc->NodeNum() = index; //store node number in reference configuration (local == global) for const int& access of virtual function nodenum in referenceframe2D

	searchtree.AddItem(Box3D(n->Pos(),n->Pos()), index);

	return index;
	
}

HRESULT CVWScale3DTool::OnDoneScaling()
{
	HRESULT hr = S_OK;
	POSITION pos;
	CTranslate3DObject * pCTrans = NULL;
	IVector			*pvPos = NULL;

	for( pos = m_TransformList.GetHeadPosition(); pos != NULL; )
	{
		pCTrans = m_TransformList.GetNext( pos );
		if(NULL != pCTrans)
		{
			pvPos = NULL;
			if (FAILED(hr = pCTrans->m_pThing->get_ObjectProperty(bstrScale, (IObjectProperty**)&pvPos)))
				goto FAIL_EXIT;

			pvPos->set(pCTrans->currentLocation.x, pCTrans->currentLocation.y, pCTrans->currentLocation.z);
		}
	}

FAIL_EXIT:
	SAFERELEASE(pvPos);
	return hr;
}

void SequenceToBatch::sequence2BatchCopy(Matrix &batch,
                                         Matrix &sequence,
                                         IVector &seq2BatchIdx,
                                         bool seq2batch) {
  int seqWidth = sequence.getWidth();
  int batchCount = batch.getHeight();
  real *batchData = batch.getData();
  real *seqData = sequence.getData();
  int *idxData = seq2BatchIdx.getData();

  if (useGpu_) {
    hl_sequence2batch_copy(
        batchData, seqData, idxData, seqWidth, batchCount, seq2batch);
  } else {
    if (seq2batch) {
#ifdef PADDLE_USE_MKLML
      const int blockMemSize = 8 * 1024;
      const int blockSize = blockMemSize / sizeof(real);
#pragma omp parallel for collapse(2)
      for (int i = 0; i < batchCount; ++i) {
        for (int j = 0; j < seqWidth; j += blockSize) {
          memcpy(batch.rowBuf(i) + j,
                 sequence.rowBuf(idxData[i]) + j,
                 (j + blockSize > seqWidth) ? (seqWidth - j) * sizeof(real)
                                            : blockMemSize);
        }
      }
#else
      for (int i = 0; i < batchCount; ++i) {
        memcpy(batch.rowBuf(i),
               sequence.rowBuf(idxData[i]),
               seqWidth * sizeof(real));
      }
#endif
    } else {
#ifdef PADDLE_USE_MKLML
#pragma omp parallel for
#endif
      for (int i = 0; i < batchCount; ++i) {
        memcpy(sequence.rowBuf(idxData[i]),
               batch.rowBuf(i),
               seqWidth * sizeof(real));
      }
    }
  }
}

void SequenceToBatch::sequence2BatchAdd(Matrix &batch,
                                        Matrix &sequence,
                                        IVector &seq2BatchIdx,
                                        bool seq2batch) {
  int seqWidth = sequence.getWidth();
  int batchCount = batch.getHeight();
  real *batchData = batch.getData();
  real *seqData = sequence.getData();
  int *idxData = seq2BatchIdx.getData();

  if (useGpu_) {
    hl_sequence2batch_add(
        batchData, seqData, idxData, seqWidth, batchCount, seq2batch);
  } else {
    for (int i = 0; i < batchCount; ++i) {
      if (seq2batch) {
        batch.subMatrix(i, 1)->add(*sequence.subMatrix(idxData[i], 1));
      } else {
        sequence.subMatrix(idxData[i], 1)->add(*batch.subMatrix(i, 1));
      }
    }
  }
}

boolean CAlgorithmClassifierBliffCFIS::classify(const IFeatureVector& rFeatureVector, float64& rf64Class, IVector& rClassificationValues)
{
	bliff::FeatureVector l_oFeatureVector(rFeatureVector.getSize(), 0);
	for(uint32 j=0; j<rFeatureVector.getSize(); j++)
	{
		l_oFeatureVector[j]=rFeatureVector[j];
	}

	FILE* l_pFile=::fopen(_ParameterFile_, "wb");
	::fwrite(m_oConfiguration.getDirectPointer(), m_oConfiguration.getSize(), 1, l_pFile);
	::fclose(l_pFile);

	itpp::Vec<double> l_vResult;
	double l_dResult;
	bliff::CFIS l_oBliffCFISClassifier;

	l_oBliffCFISClassifier.readParams(_ParameterFile_);
	l_vResult=l_oBliffCFISClassifier.classify(l_oFeatureVector);
	l_dResult=l_oBliffCFISClassifier.assign(l_oFeatureVector);

	if(ip_ui64OutputMode == OVP_TypeId_CFISOutputMode_ClassMembership.toUInteger())
	{
		itpp::Vec<double> l_vTmpVec(l_oBliffCFISClassifier.getNbClasses());
		l_vTmpVec.zeros();
		std::vector<double> l_oClassLabels = l_oBliffCFISClassifier.getClassLabels();

		//finding the maximal degree of fulfillment for each class
		for(uint32 j=0; j < l_oBliffCFISClassifier.getNbRules(); j++)
		{
			for(uint32 i=0; i < l_oClassLabels.size(); i++)
			{
				if((l_oBliffCFISClassifier.getRule(j)->getClass() == l_oClassLabels[i]) && (l_vResult[j]>l_vTmpVec[i]))
				{
					l_vTmpVec[i]=l_vResult[j];
				}
			}
		}

		//switching from rule fulfillment to class membership
		l_vResult=l_vTmpVec;
	}

	//converting from BLiFF++ output to OpenViBE output (classification state)
	rf64Class=l_dResult;
	rClassificationValues.setSize(l_vResult.size());
	for(size_t i=0; i<rClassificationValues.getSize(); i++)
	{
		rClassificationValues[i]=l_vResult[i];
	}

	//labelling the element of the output vector
	if(ip_ui64OutputMode == OVP_TypeId_CFISOutputMode_ClassMembership.toUInteger())
	{
		char l_sBuffer[1024];
		std::vector<double> l_vClassLabels = l_oBliffCFISClassifier.getClassLabels();
		for(uint32 i=0; i < rClassificationValues.getSize(); i++)
		{
			sprintf(l_sBuffer, "Class %d membership degree", (uint32)l_vClassLabels[i]);
			rClassificationValues.setElementLabel(i, l_sBuffer);
		}
	}
	else if(ip_ui64OutputMode == OVP_TypeId_CFISOutputMode_RuleFulfillment.toUInteger())
	{
		char l_sBuffer[1024];
		for(uint32 i=0; i < rClassificationValues.getSize(); i++)
		{
			sprintf(l_sBuffer, "Rule %d degree of fulfillment", (uint32)i+1);
			rClassificationValues.setElementLabel(i, l_sBuffer);
		}
	}
	else
	{
		this->getLogManager() << LogLevel_Warning << "Unhandled CFIS output mode " << ip_ui64OutputMode << " (" << this->getTypeManager().getEnumerationEntryNameFromValue(OVP_TypeId_CFISOutputMode, ip_ui64OutputMode) << "\n";
	}

	return true;
}