| 44 | 
  | 
 * | 
| 45 | 
  | 
 *  Created by Charles F. Vardeman II on 11 Dec 2006. | 
| 46 | 
  | 
 *  @author  Charles F. Vardeman II | 
| 47 | 
< | 
 *  @version $Id: ConvexHull.cpp,v 1.10 2008-10-15 18:26:01 chuckv Exp $ | 
| 47 | 
> | 
 *  @version $Id: ConvexHull.cpp,v 1.12 2008-10-21 16:44:00 chuckv Exp $ | 
| 48 | 
  | 
 * | 
| 49 | 
  | 
 */ | 
| 50 | 
  | 
 | 
| 370 | 
  | 
#ifdef HAVE_QHULL | 
| 371 | 
  | 
/* Old options Qt Qu Qg QG0 FA */ | 
| 372 | 
  | 
/* More old opts Qc Qi Pp*/ | 
| 373 | 
< | 
ConvexHull::ConvexHull() : Hull(), dim_(3), options_("qhull Qt Pp"), Ns_(200) { | 
| 373 | 
> | 
ConvexHull::ConvexHull() : Hull(), dim_(3), options_("qhull Qt Pp"), Ns_(200), nTriangles_(0) { | 
| 374 | 
  | 
  //If we are doing the mpi version, set up some vectors for data communication | 
| 375 | 
  | 
#ifdef IS_MPI | 
| 376 | 
  | 
 | 
| 438 | 
  | 
   | 
| 439 | 
  | 
  boolT ismalloc = False; | 
| 440 | 
  | 
  /* Clean up memory from previous convex hull calculations*/ | 
| 441 | 
+ | 
   | 
| 442 | 
  | 
  Triangles_.clear(); | 
| 443 | 
  | 
  surfaceSDs_.clear(); | 
| 444 | 
  | 
  surfaceSDs_.reserve(Ns_); | 
| 447 | 
  | 
                    const_cast<char *>(options_.c_str()), NULL, stderr)) { | 
| 448 | 
  | 
 | 
| 449 | 
  | 
      sprintf(painCave.errMsg, "ConvexHull: Qhull failed to compute convex hull"); | 
| 450 | 
< | 
      painCave.isFatal = 0; | 
| 450 | 
> | 
      painCave.isFatal = 1; | 
| 451 | 
  | 
      simError(); | 
| 452 | 
  | 
       | 
| 453 | 
  | 
  } //qh_new_qhull | 
| 455 | 
  | 
 | 
| 456 | 
  | 
#ifdef IS_MPI | 
| 457 | 
  | 
  std::vector<double> localPts; | 
| 458 | 
+ | 
  std::vector<double> localVel; | 
| 459 | 
  | 
  int localPtArraySize; | 
| 460 | 
  | 
   | 
| 461 | 
  | 
  | 
| 467 | 
  | 
    if (!facet->simplicial){ | 
| 468 | 
  | 
      // should never happen with Qt | 
| 469 | 
  | 
      sprintf(painCave.errMsg, "ConvexHull: non-simplicaial facet detected"); | 
| 470 | 
< | 
      painCave.isFatal = 0; | 
| 470 | 
> | 
      painCave.isFatal = 1; | 
| 471 | 
  | 
      simError(); | 
| 472 | 
  | 
    } | 
| 473 | 
  | 
     | 
| 519 | 
  | 
    localPts.push_back(ptArray[dim_ * idx]);      | 
| 520 | 
  | 
    localPts.push_back(ptArray[dim_ * idx + 1]);  | 
| 521 | 
  | 
    localPts.push_back(ptArray[dim_ * idx + 2]); | 
| 522 | 
+ | 
 | 
| 523 | 
+ | 
    Vector3d vel = bodydoubles[idx]->getVel(); | 
| 524 | 
+ | 
    localVel.push_back(vel.x()); | 
| 525 | 
+ | 
    localVel.push_back(vel.y()); | 
| 526 | 
+ | 
    localVel.push_back(vel.z()); | 
| 527 | 
+ | 
 | 
| 528 | 
  | 
    localPtsMap.push_back(idx);  | 
| 529 | 
  | 
  } | 
| 530 | 
  | 
 | 
| 543 | 
  | 
  int nglobalPts = int(Nsglobal_/3); | 
| 544 | 
  | 
  | 
| 545 | 
  | 
 | 
| 546 | 
< | 
  std::vector<double> globalPts; | 
| 547 | 
< | 
  globalPts.resize(Nsglobal_); | 
| 546 | 
> | 
  std::vector<double> globalPts(Nsglobal_); | 
| 547 | 
> | 
  std::vector<double> globalVel(Nsglobal_); | 
| 548 | 
  | 
 | 
| 549 | 
  | 
  isSurfaceID.resize(nglobalPts); | 
| 550 | 
  | 
 | 
| 562 | 
  | 
  /* gather the potential hull */ | 
| 563 | 
  | 
   | 
| 564 | 
  | 
  MPI::COMM_WORLD.Allgatherv(&localPts[0],localPtArraySize,MPI::DOUBLE,&globalPts[0],&NstoProc_[0],&displs_[0],MPI::DOUBLE); | 
| 565 | 
+ | 
  MPI::COMM_WORLD.Allgatherv(&localVel[0],localPtArraySize,MPI::DOUBLE,&globalVel[0],&NstoProc_[0],&displs_[0],MPI::DOUBLE); | 
| 566 | 
  | 
 | 
| 567 | 
  | 
  /* | 
| 568 | 
  | 
  if (myrank_ == 0){ | 
| 609 | 
  | 
        simError(); | 
| 610 | 
  | 
      } //simplicical | 
| 611 | 
  | 
       | 
| 612 | 
< | 
      Triangle* face = new Triangle(); | 
| 612 | 
> | 
      Triangle face; | 
| 613 | 
  | 
      Vector3d  V3dNormal(facet->normal[0],facet->normal[1],facet->normal[2]); | 
| 614 | 
< | 
      face->setNormal(V3dNormal); | 
| 614 | 
> | 
      face.setNormal(V3dNormal); | 
| 615 | 
  | 
  | 
| 616 | 
  | 
       | 
| 617 | 
  | 
 | 
| 618 | 
  | 
      RealType faceArea = 0.5*V3dNormal.length(); | 
| 619 | 
< | 
      face->setArea(faceArea); | 
| 619 | 
> | 
      face.setArea(faceArea); | 
| 620 | 
  | 
 | 
| 621 | 
  | 
 | 
| 622 | 
  | 
      vertices = qh_facet3vertex(facet); | 
| 623 | 
  | 
       | 
| 624 | 
  | 
      coordT *center = qh_getcenter(vertices); | 
| 625 | 
  | 
      Vector3d V3dCentroid(center[0], center[1], center[2]); | 
| 626 | 
< | 
      face->setCentroid(V3dCentroid); | 
| 627 | 
< | 
 | 
| 626 | 
> | 
      face.setCentroid(V3dCentroid); | 
| 627 | 
> | 
      Vector3d faceVel = V3Zero; | 
| 628 | 
  | 
      FOREACHvertex_(vertices){ | 
| 629 | 
  | 
        id = qh_pointid(vertex->point); | 
| 630 | 
  | 
        int localindex = id; | 
| 631 | 
  | 
#ifdef IS_MPI | 
| 632 | 
< | 
         | 
| 632 | 
> | 
        Vector3d velVector(globalVel[dim_ * id],globalVel[dim_ * id + 1], globalVel[dim_ * id + 1]); | 
| 633 | 
> | 
        faceVel = faceVel + velVector; | 
| 634 | 
  | 
        if (id >= noffset/3 && id < (noffset + localPtArraySize)/3 ){ | 
| 635 | 
  | 
          localindex = localPtsMap[id-noffset/3]; | 
| 636 | 
+ | 
#else | 
| 637 | 
+ | 
          faceVel = faceVel + bodydoubles[localindex]->getVel(); | 
| 638 | 
  | 
#endif | 
| 639 | 
< | 
          face->addVertex(bodydoubles[localindex]); | 
| 639 | 
> | 
          face.addVertex(bodydoubles[localindex]); | 
| 640 | 
  | 
          if( !isSurfaceID[id] ){ | 
| 641 | 
  | 
            isSurfaceID[id] = true; | 
| 642 | 
  | 
#ifdef IS_MPI        | 
| 650 | 
  | 
#ifdef IS_MPI | 
| 651 | 
  | 
          | 
| 652 | 
  | 
        }else{ | 
| 653 | 
< | 
          face->addVertex(NULL); | 
| 653 | 
> | 
          face.addVertex(NULL); | 
| 654 | 
  | 
          } | 
| 655 | 
  | 
#endif | 
| 656 | 
  | 
      } //Foreachvertex | 
| 661 | 
  | 
          surfaceSDs_.push_back(bodydoubles[id]); | 
| 662 | 
  | 
        } | 
| 663 | 
  | 
      } | 
| 664 | 
< | 
 | 
| 664 | 
> | 
      */ | 
| 665 | 
> | 
      face.setFacetVelocity(faceVel/3.0); | 
| 666 | 
  | 
      Triangles_.push_back(face); | 
| 667 | 
  | 
      qh_settempfree(&vertices);       | 
| 668 | 
< | 
      */ | 
| 668 | 
> | 
 | 
| 669 | 
  | 
    } //FORALLfacets | 
| 670 | 
  | 
 | 
| 671 | 
  | 
    /* | 
| 679 | 
  | 
 | 
| 680 | 
  | 
 | 
| 681 | 
  | 
    Ns_ = surfaceSDs_.size(); | 
| 682 | 
+ | 
    nTriangles_ = Triangles_.size(); | 
| 683 | 
  | 
     | 
| 670 | 
– | 
     | 
| 684 | 
  | 
    qh_getarea(qh facet_list); | 
| 685 | 
  | 
    volume_ = qh totvol; | 
| 686 | 
  | 
    area_ = qh totarea; |