directions3D.cpp

// Gmsh - Copyright (C) 1997-2013 C. Geuzaine, J.-F. Remacle
//
// See the LICENSE.txt file for license information. Please report all
// bugs and problems to the public mailing list <gmsh@geuz.org>.
//
// Contributor(s):
//   Tristan Carrier François Henrotte

#include <fstream>
#include "GModel.h"
#include "BackgroundMesh.h"
#include "meshGFaceDelaunayInsertion.h"
#include "MTetrahedron.h"
#include "directions3D.h"

#if defined(HAVE_PETSC)
#include "dofManager.h"
#include "laplaceTerm.h"
#include "linearSystemPETSc.h"
#include "linearSystemFull.h"
#endif



/****************class Frame_field****************/

Frame_field::Frame_field(){}

void Frame_field::init_region(GRegion* gr){
  // Fill in a ANN tree with the bondary cross field of region gr
#if defined(HAVE_ANN)
  int index;
  MVertex* vertex;
  GFace* gf;
  std::list<GFace*> faces;
  STensor3 m(1.0);

  Nearest_point::init_region(gr);

  faces = gr->faces();

  temp.clear();
  field.clear();

  for( std::list<GFace*>::iterator it=faces.begin();it!=faces.end();it++){
    gf = *it;
    init_face(gf);
  }

  ANNpointArray duplicate = annAllocPts(temp.size(),3);

  index = 0;
  for(std::map<MVertex*,STensor3>::iterator it=temp.begin(); it != temp.end(); it++){
    vertex = it->first;
    m = it->second;
    duplicate[index][0] = vertex->x();
    duplicate[index][1] = vertex->y();
    duplicate[index][2] = vertex->z();
    field.push_back(std::pair<MVertex*,STensor3>(vertex,m));
    index++;
  }

  kd_tree = new ANNkd_tree(duplicate, temp.size(), 3);
#endif
}

void Frame_field::init_face(GFace* gf){
  // Fills the auxiliary std::map "temp" with a pair <vertex, STensor3>
  // for each vertex of the face gf.
  unsigned int i;