C++ CollOfScalar类代码示例

CollOfScalar EquelleRuntimeCUDA::operatorOn(const CollOfScalar& data_in,
					    const CollOfIndices<codim>& from_set,
					    const CollOfIndices<codim>& to_set) {

    if ( data_in.size() != from_set.size()) {
	OPM_THROW(std::logic_error, "data_in (size " << data_in.size() << ") and from_set (size " << from_set.size() << ") have to be of the same size in On function for CollOfScalar.");
    }
    //if ( to_set.size() > from_set.size() ) {
    //	OPM_THROW(std::runtime_error, "To_set (size " << to_set.size() << ") has to be a subset of from_set (size " << from_set.size() << ")");
    //}

    return dev_grid_.operatorOn(data_in, from_set, to_set);
}

CollOfScalar EquelleRuntimeCUDA::operatorExtend(const CollOfScalar& data_in,
						const CollOfIndices<codim>& from_set,
						const CollOfIndices<codim>& to_set) const {
    
    if (data_in.size() != from_set.size() ) {
	OPM_THROW(std::runtime_error, "data_in (size " << data_in.size() << ") and from_set (size " << from_set.size() << ") have to be of the same size in Extend function."); 
    }
    if (from_set.size() > to_set.size() ) {
	OPM_THROW(std::runtime_error, "From_set (size " << from_set.size() << ") has to be a subset of to_set (size " << to_set.size() << ")");
    }
    
    return dev_grid_.operatorExtend(data_in, from_set, to_set);
}

int compare(CollOfScalar scal, double sol[], 
	    int sol_size,
	    std::string test) 
{ 
    // Test size:
    if ( scal.size() != sol_size ) {
	std::cout << "Error in valsOnGrid.cpp - testing " << test << "\n";
	std::cout << "\tThe collection is of wrong size!\n";
	std::cout << "\tSize is " << scal.size() << " but should be " << sol_size << "\n";
	return 1;
    }
    
    // Testing indices
    std::vector<double> host = scal.copyToHost();
    std::cout << "CollOfScalar " << test << " is the following:\n";
    bool correct = true;
    for (int i = 0; i < host.size(); ++i) {
	std::cout << host[i] << " ";
	if (i < sol_size) {
	    //if (host[i] != sol[i]) {
	    if ( fabs(host[i] - sol[i]) > 1000*std::numeric_limits<double>::epsilon() ) {
		std::cout << "(<- " << sol[i] << ") ";
		correct = false;
	    }
	}
    }
    if (correct) {
	std::cout << "\n\tThis is correct\n";
    } else {
	std::cout << "\n\tThis is wrong\n";
	std::cout << "Error in valsOnGrid.cpp - testing " << test << "\n";
	std::cout << "\tThe indices in the collection is wrong\n";
	return 1;
    }

    return 0;

}

CollOfScalar EquelleRuntimeCUDA::newtonSolve(const ResidualFunctor& rescomp,
                                            const CollOfScalar& u_initialguess)
{
    Opm::time::StopWatch clock;
    clock.start();

    // Set up Newton loop.
 
    // Define the primary variable
    CollOfScalar u = CollOfScalar(u_initialguess, true);
 
    if (verbose_ > 2) {
        output("Initial u", u);
        output("    newtonSolve: norm (initial u)", twoNorm(u));
    }
    CollOfScalar residual = rescomp(u);   
    if (verbose_ > 2) {
        output("Initial residual", residual);
        output("    newtonSolve: norm (initial residual)", twoNorm(residual));
    }

    int iter = 0;

    // Debugging output not specified in Equelle.
    if (verbose_ > 1) {
        std::cout << "    newtonSolve: iter = " << iter << " (max = " << max_iter_
		  << "), norm(residual) = " << twoNorm(residual)
                  << " (tol = " << abs_res_tol_ << ")" << std::endl;
    }

    CollOfScalar du;

    // Execute newton loop until residual is small or we have used too many iterations.
    while ( (twoNorm(residual) > abs_res_tol_) && (iter < max_iter_) ) {
	
	if ( solver_.getSolver() == CPU ) {
	    du = serialSolveForUpdate(residual);
	}
	else {
	    // Solve linear equations for du, apply update.
	    du = solver_.solve(residual.derivative(),
			       residual.value(),
			       verbose_);
	}

	// du is a constant, hence, u is still a primary variable with an identity
	// matrix as its derivative.
	u = u - du;

        // Recompute residual.
        residual = rescomp(u);

        if (verbose_ > 2) {
            // Debugging output not specified in Equelle.
            output("u", u);
            output("    newtonSolve: norm(u)", twoNorm(u));
            output("residual", residual);
            output("    newtonSolve: norm(residual)", twoNorm(residual));
        }

        ++iter;

        // Debugging output not specified in Equelle.
        if (verbose_ > 1) {
            std::cout << "    newtonSolve: iter = " << iter << " (max = " << max_iter_
		      << "), norm(residual) = " << twoNorm(residual)
                      << " (tol = " << abs_res_tol_ << ")" << std::endl;
        }

    }
    if (verbose_ > 0) {
        if (twoNorm(residual) > abs_res_tol_) {
            std::cout << "Newton solver failed to converge in " << max_iter_ << " iterations" << std::endl;
        } else {
            std::cout << "Newton solver converged in " << iter << " iterations" << std::endl;
        }
    }

    if (verbose_ > 1) {
        std::cout << "Newton solver took: " << clock.secsSinceLast() << " seconds." << std::endl;
    }

    return CollOfScalar(u.value());
}