C++ Param::boundingBox方法代码示例

// Write VTK file for the flow around the obstacle, to be viewed with Paraview.
void writeVTK(OffLatticeBoundaryCondition3D<T,DESCRIPTOR,Velocity>& bc, plint iT)
{
    VtkImageOutput3D<T> vtkOut(createFileName("volume", iT, PADDING));
    vtkOut.writeData<float>( *bc.computeVelocityNorm(param.boundingBox()),
                             "velocityNorm", param.dx/param.dt );
    vtkOut.writeData<3,float>(*bc.computeVelocity(param.boundingBox()), "velocity", param.dx/param.dt);
    vtkOut.writeData<float>( *bc.computePressure(param.boundingBox()),
                             "pressure", param.dx*param.dx/(param.dt*param.dt) );
}

void runProgram()
{
    /*
     * Read the obstacle geometry.
     */

    pcout << std::endl << "Reading STL data for the obstacle geometry." << std::endl;
    Array<T,3> center(param.cx, param.cy, param.cz);
    Array<T,3> centerLB(param.cxLB, param.cyLB, param.czLB);
    // The triangle-set defines the surface of the geometry.
    TriangleSet<T> triangleSet(param.geometry_fname, DBL);

    // Place the obstacle in the correct place in the simulation domain.
    // Here the "geometric center" of the obstacle is computed manually,
    // by computing first its bounding cuboid. In cases that the STL
    // file with the geometry of the obstacle contains its center as
    // the point, say (0, 0, 0), then the following variable
    // "obstacleCenter" must be set to (0, 0, 0) manually.
    Cuboid<T> bCuboid = triangleSet.getBoundingCuboid();
    Array<T,3> obstacleCenter = 0.5 * (bCuboid.lowerLeftCorner + bCuboid.upperRightCorner);
    triangleSet.translate(-obstacleCenter);
    triangleSet.scale(1.0/param.dx); // In lattice units from now on...
    triangleSet.translate(centerLB);
    triangleSet.writeBinarySTL(outputDir+"obstacle_LB.stl");

    // The DEFscaledMesh, and the triangle-boundary are more sophisticated data
    // structures used internally by Palabos to treat the boundary.
    plint xDirection = 0;
    plint borderWidth = 1;      // Because Guo acts in a one-cell layer.
                                // Requirement: margin>=borderWidth.
    plint margin = 1;           // Extra margin of allocated cells around the obstacle, for the case of moving walls.
    plint blockSize = 0;        // Size of blocks in the sparse/parallel representation.
                                // Zero means: don't use sparse representation.
    DEFscaledMesh<T> defMesh(triangleSet, 0, xDirection, margin, Dot3D(0, 0, 0));
    TriangleBoundary3D<T> boundary(defMesh);
    //boundary.getMesh().inflate();

    pcout << "tau = " << 1.0/param.omega << std::endl;
    pcout << "dx = " << param.dx << std::endl;
    pcout << "dt = " << param.dt << std::endl;
    pcout << "Number of iterations in an integral time scale: " << (plint) (1.0/param.dt) << std::endl;

    /*
     * Voxelize the domain.
     */

    // Voxelize the domain means: decide which lattice nodes are inside the obstacle and which are outside.
    pcout << std::endl << "Voxelizing the domain." << std::endl;
    plint extendedEnvelopeWidth = 2;   // Extrapolated off-lattice BCs.
    const int flowType = voxelFlag::outside;
    VoxelizedDomain3D<T> voxelizedDomain (
            boundary, flowType, param.boundingBox(), borderWidth, extendedEnvelopeWidth, blockSize );
    pcout << getMultiBlockInfo(voxelizedDomain.getVoxelMatrix()) << std::endl;

    /*
     * Generate the lattice, the density and momentum blocks.
     */

    pcout << "Generating the lattice, the rhoBar and j fields." << std::endl;
    MultiBlockLattice3D<T,DESCRIPTOR> *lattice = new MultiBlockLattice3D<T,DESCRIPTOR>(voxelizedDomain.getVoxelMatrix());
    if (param.useSmago) {
        defineDynamics(*lattice, lattice->getBoundingBox(),
                new SmagorinskyBGKdynamics<T,DESCRIPTOR>(param.omega, param.cSmago));
        pcout << "Using Smagorinsky BGK dynamics." << std::endl;
    } else {
        defineDynamics(*lattice, lattice->getBoundingBox(),
                new BGKdynamics<T,DESCRIPTOR>(param.omega));
        pcout << "Using BGK dynamics." << std::endl;
    }
    bool velIsJ = false;
    defineDynamics(*lattice, voxelizedDomain.getVoxelMatrix(), lattice->getBoundingBox(),
            new NoDynamics<T,DESCRIPTOR>(), voxelFlag::inside);
    lattice->toggleInternalStatistics(false);

    MultiBlockManagement3D sparseBlockManagement(lattice->getMultiBlockManagement());

    // The rhoBar and j fields are used at both the collision and at the implementation of the
    // outflow boundary condition.
    plint envelopeWidth = 1;
    MultiScalarField3D<T> *rhoBar = new MultiScalarField3D<T> (
            MultiBlockManagement3D (
                sparseBlockManagement.getSparseBlockStructure(),
                sparseBlockManagement.getThreadAttribution().clone(),
                envelopeWidth ),
            defaultMultiBlockPolicy3D().getBlockCommunicator(),
            defaultMultiBlockPolicy3D().getCombinedStatistics(),
            defaultMultiBlockPolicy3D().getMultiScalarAccess<T>() );
    rhoBar->toggleInternalStatistics(false);

    MultiTensorField3D<T,3> *j = new MultiTensorField3D<T,3> (
            MultiBlockManagement3D (
                sparseBlockManagement.getSparseBlockStructure(),
                sparseBlockManagement.getThreadAttribution().clone(),
                envelopeWidth ),
            defaultMultiBlockPolicy3D().getBlockCommunicator(),
            defaultMultiBlockPolicy3D().getCombinedStatistics(),
            defaultMultiBlockPolicy3D().getMultiTensorAccess<T,3>() );
    j->toggleInternalStatistics(false);

    std::vector<MultiBlock3D*> lattice_rho_bar_j_arg;
//.........这里部分代码省略.........