C++ GridFunction::GetGridFunction方法代码示例

int IO::getAmountOfFluidcells(GridFunction& geo)
{
	int aof = 0;
	for (int i = geo.beginwrite[0]; i <= geo.endwrite[0]; i++)
	{
		for  (int j = geo.beginwrite[1]; j <= geo.endwrite[1]; j++ )
		{
			if(geo.GetGridFunction(i,j) >= 16)
			    aof++;
		}
	}
	return aof;
}

RealType Solver::computeResidual(GridFunction& sourcegridfunction,
    				     GridFunctionType& rhs,
    					 const PointType& h){
//<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
    // The pre-value to be returned (return sqrt(doubleSum)):
    RealType doubleSum = 0.0;

    /* We need to compute the derivatives p_xx and p_yy, therefore the stencil has to be applied.
     */

    MultiIndexType dim = sourcegridfunction.GetGridDimension();
    MultiIndexType bread  (0,0);
    MultiIndexType eread  (dim[0]-1,dim[1]-1);
    MultiIndexType bwrite (1,1);
    MultiIndexType ewrite (dim[0]-2,dim[1]-2);

    //Compute the needed derivations for the whole (inner?) area
    Stencil stencil(3,h); 					// bzw. Kann man einfach const weitergeben? /Wie?
    //Get the values for derivative in x-direction:
    GridFunction Fxx(dim);
    stencil.ApplyFxxStencilOperator(bread, eread, bwrite, ewrite, sourcegridfunction.GetGridFunction(), Fxx);
    //Get the values for derivative in y-direction:
    GridFunction Fyy(dim);
    stencil.ApplyFyyStencilOperator(bread, eread, bwrite, ewrite, sourcegridfunction.GetGridFunction(), Fyy);

    // Compute the residual: res = sqrt(Sum_i^I(Sum_j^J((p_xx+p_yy-rightHandSide)²/(I*J))))
    RealType derivator;
    for (IndexType i = 1; i <= dim[0]-2; i++)
    {
    	for (IndexType j = 1; j <= dim[1]-2; j++)
    	{
    		derivator = Fxx.GetGridFunction()[i][j]+ Fyy.GetGridFunction()[i][j] - rhs[i][j];
            doubleSum +=  derivator*derivator / (dim[0]-2) / (dim[1]-2);
    	}
    }
    //std::cout<<doubleSum<<std::endl;
    return sqrt(doubleSum);
}

void IO::writeVTKSlavefile(GridFunction& u_gridfunction,
		  GridFunction& v_gridfunction, GridFunction& p_gridfunction,
		  GridFunction& T_gridfunction, GridFunction& Geo_gridfunction,
		  const PointType& delta, int mpiSizeH, int mpiSizeV, int step,
		  int rank){

	double deltaX =delta[0];
	double deltaY =delta[1];
	int ibegin = p_gridfunction.beginwrite[0];
	int iend   = p_gridfunction.endwrite[0];
	int jbegin = p_gridfunction.beginwrite[1];
	int jend   = p_gridfunction.endwrite[1];
	GridFunctionType p = p_gridfunction.GetGridFunction();
	GridFunctionType u = u_gridfunction.GetGridFunction();
	GridFunctionType v = v_gridfunction.GetGridFunction();
	GridFunctionType T = T_gridfunction.GetGridFunction();
	GridFunctionType geo = Geo_gridfunction.GetGridFunction();
	int localgriddimensionX = iend-ibegin+1;
	int localgriddimensionY = jend-jbegin+1;

	int processorgridcoordX = rank % mpiSizeH;
	int processorgridcoordY = floor(rank / mpiSizeH);

	int x1=processorgridcoordX    *localgriddimensionX-processorgridcoordX;
	int x2=(processorgridcoordX+1)*localgriddimensionX-processorgridcoordX-1;
	int x3=processorgridcoordY    *localgriddimensionY-processorgridcoordY;
	int x4=(processorgridcoordY+1)*localgriddimensionY-processorgridcoordY-1;

	  char numstr[21];
	  sprintf (numstr, "%d", step);
	  std::string filename;
	  filename.append ("./");
	  filename.append (output);
	  filename.append ("/");
	  filename.append ("sol_");
	  filename.append (numstr);
	  filename.append ("_rank");
	  sprintf (numstr, "%d", rank);
	  filename.append (numstr);
	  filename.append (".vtr");

	  std::filebuf fb;
	  fb.open (const_cast < char *>(filename.c_str ()), std::ios::out);
	  std::ostream os (&fb);
	  os << "<?xml version=\"1.0\"?>" << std::endl
	  //  << "<VTKFile type=\"RectilinearGrid\">" << std::endl
		<< "<VTKFile type=\"RectilinearGrid\">" << std::endl
	    << "<RectilinearGrid WholeExtent=\""
	    << x1 << " " << x2 << " "
	    << x3 << " " << x4 << " "
	    << "0" << " " << "0" << " "
	    << "\" GhostLevel=\"" << "0" << "\">" << std::endl
	    << "<Piece Extent=\""<<x1<<" "<<x2<<" "<<x3<<" "<<x4<<" 0 0 \">" <<std::endl
	    << "<Coordinates>"<<std::endl
	        <<
	        "<DataArray type=\"Float64\" format=\"ascii\"> "
	        << std::endl;
	      //coordinates
	      for (int i = ibegin; i <= iend; ++i)
	      {
	    	  os << std::scientific <<(processorgridcoordX*localgriddimensionX)*deltaX-processorgridcoordX*deltaX+ i * deltaX << " ";
	      }
	      os << std::endl;
	      os << "</DataArray>" << std::endl
	      << "<DataArray type=\"Float64\" format=\"ascii\"> " << std::endl;
	      for (int j = jbegin; j <= jend; ++j)
	      {
	    	  os << std::scientific << (processorgridcoordY*localgriddimensionY)*deltaY-processorgridcoordY*deltaY + j * deltaY << " ";
	      }
	      os << std::endl
	      << "</DataArray>" << std::endl
	      << "<DataArray type=\"Float64\" format=\"ascii\"> " << std::endl
	      << "0 0" << std::endl
	      << "</DataArray>" << std::endl

	        << "</Coordinates>" << std::endl
	        << "<PointData>"
	        << std::endl <<
	        "<DataArray Name=\"field\" NumberOfComponents=\"3\" type=\"Float64\" >" <<
	        std::endl;
	      //velocities
	      for (int j = jbegin; j <= jend; ++j)
	      	    	{
	      	    	  RealType y = j * deltaY;

	      for (int i = ibegin; i <= iend; ++i)
	        {
	          RealType x = i * deltaX;



	    	  os << std::scientific << interpolateVelocityU (x, y, u,
	    							 delta) << " " <<
	    	    interpolateVelocityV (x, y, v, delta) << " " << 0. << std::endl;
	    	}

	        }
	      os << "</DataArray>" << std::endl
	        << "<DataArray type=\"Float64\" Name=\"P\" format=\"ascii\">" <<
	        std::endl;
//.........这里部分代码省略.........