C++ MyVector::Size方法代码示例

int main() {

    const std::string Help =
        "-------------------------------------------------------------------------\n"
        "TestComputeAKV:                                                          \n"
        "-------------------------------------------------------------------------\n"
        "OPTIONS:                                                                 \n"
        "Nth=<int>  theta resolution [default 3]                                  \n"
        "Nph=<int>  phi resolution [default 4]                                    \n"
        "Radius=<double> radius of sphere. [default 1.0]                          \n"
        "AKVGuess=MyVector<double> a guess for the values of THETA, thetap, phip  \n"
        "         [default (0.,0.,0.)]                                            \n"
        "L_resid_tol=<double> tolerance for L residuals when finding approximate  \n"
        "            Killing vectors.  [default 1.e-12]                           \n"
        "v_resid_tol=<double> tolerance for v residuals when finding approximate  \n"
        "            Killing vectors.  [default 1.e-12]                           \n"
        "min_thetap = for values less than this, thetap is considered close to    \n"
        "               zero. [default 1.e-5]                                     \n"
        "symmetry_tol=<double> abs(THETA) must be less than this value to be      \n"
        "             considered an exact symmetry.  [default 1.e-11]             \n"
        "ResidualSize=<double> determines the tolerance for residuals from the    \n"
        "             multidimensional root finder.  [default to 1.e-11]          \n"
        "Solver = <std::string> which gsl multidimensional root finding algorith  \n"
        "        should be used. [default Newton]                                 \n"
        "Verbose=<bool> Print spectral coefficients and warnings if true          \n"
        "        [default false]                                                  \n"
        ;

    std::string Options = ReadFileIntoString("Test.input");
    OptionParser op(Options,Help);
    const int Nth = op.Get<int>("Nth", 3);
    const int Nph = op.Get<int>("Nph", 4);
    const double rad = op.Get<double>("Radius",1.0);
    MyVector<double> AKVGuess =
        op.Get<MyVector<double> >("AKVGuess",MyVector<double>(MV::Size(3),0.0));
    //must be three-dimensional
    REQUIRE(AKVGuess.Size()==3,"AKVGuess has Size " << AKVGuess.Size()
            << ", should be 3.");
    const double L_resid_tol = op.Get<double>("L_resid_tol", 1.e-12);
    const double v_resid_tol = op.Get<double>("L_resid_tol", 1.e-12);
    const double residualSize = op.Get<double>("ResidualSize", 1.e-11);
    const double min_thetap = op.Get<double>("min_theta",1.e-5);
    const double symmetry_tol = op.Get<double>("symmetry_tol",1.e-11);
    const std::string solver = op.Get<std::string>("Solver","Newton");
    const bool verbose = op.Get<bool>("Verbose", false);
    const MyVector<bool> printDiagnostic = MyVector<bool>(MV::Size(6), true);

    //create skm
    const StrahlkorperMesh skm(Nth, Nph);
    //create surface basis
    const SurfaceBasis sb(skm);
    //get theta, phi
    const DataMesh theta(skm.SurfaceCoords()(0));
    const DataMesh phi(skm.SurfaceCoords()(1));

    //set the initial guesses to be along particular axes
    const int axes = 3; //the number of perpendicular axes

    //create conformal factors for every rotation
    const int syms = 5; //the number of axisymmetries we are testing

    for(int s=4; s<5; s++) { //index over conformal factor symmetries
        //for(int s=0; s<syms; s++){//index over conformal factor symmetries
        //create conformal factor
        const DataMesh Psi = ConstructConformalFactor(theta, phi, s);

        //set the initial guesses
        double THETA[3] = {AKVGuess[0],0.,0.};
        double thetap[3] = {AKVGuess[1],0.,0.};
        double phip[3] = {AKVGuess[2],0.,0.};

        //save the v, xi solutions along particular axes
        MyVector<DataMesh> v(MV::Size(3),DataMesh::Empty);
        MyVector<DataMesh> rotated_v(MV::Size(3),DataMesh::Empty);
        MyVector<Tensor<DataMesh> > xi(MV::Size(axes),Tensor<DataMesh>(2,"1",DataMesh::Empty));

        //save the <v_i|v_j> inner product solutions
        double v0v0 = 0.;
        double v1v1 = 0.;
        double v2v2 = 0.;
        double v0v1 = 0.;
        double v0v2 = 0.;
        double v1v2 = 0.;
        //int symmetries[3] = 0; //counts the number of symmetries

        //compute some useful quantities
        const DataMesh rp2 = rad * Psi * Psi;
        const DataMesh r2p4 = rp2*rp2;
        const DataMesh llncf = sb.ScalarLaplacian(log(Psi));
        const DataMesh Ricci = 2.0 * (1.0-2.0*llncf) / r2p4;
        const Tensor<DataMesh> GradRicci = sb.Gradient(Ricci);

        for(int a=0; a<axes; a++) { //index over perpendicular axes to find AKV solutions

            //if the diagnostics below decide that there is a bad solution for v[a]
            //(usually a repeated solution), this flag will indicate that the
            //solver should be run again
            bool badAKVSolution = false;

            //generate a guess for the next axis of symmetry based on prior solutions.
//.........这里部分代码省略.........