Commit 7fa34b10 authored by Kirill Terekhov's avatar Kirill Terekhov
Browse files

Merge branch 'master' of https://github.com/INMOST-DEV/INMOST

parents 5f2777f8 603c8dd8
......@@ -47,8 +47,34 @@ int main(int argc, char ** argv)
std::string fout = "grid_out.pmf";
int is2D = 1;
double skew = 0;
if( argc > 2 ) is2D = atoi(argv[2]);
if( argc > 3 ) fout = std::string(argv[3]);
if( argc > 4 ) skew = atof(argv[4]);
if( skew > 0)
{
double cmax[3] = {-1.0e20,-1.0e20,-1.0e20}, cmin[3] = {1.0e20,1.0e20,1.0e20};
for(Mesh::iteratorNode n = m->BeginNode(); n != m->EndNode(); ++n)
{
Storage::real_array c = n->Coords();
for(int k = 0; k < 3; ++k)
{
if( cmax[k] < c[k] ) cmax[k] = c[k];
if( cmin[k] > c[k] ) cmin[k] = c[k];
}
}
for(Mesh::iteratorNode n = m->BeginNode(); n != m->EndNode(); ++n)
{
Storage::real_array c = n->Coords();
for(int k = 0; k < 3; ++k)
{
double t = (c[k]-cmin[k])/(cmax[k]-cmin[k]);
double a = 0.5 + 0.5*(t-0.5)/sqrt(pow(t-0.5,2)+pow(skew,2))*sqrt(0.25+pow(skew,2))*2;
c[k] = cmin[k] + (cmax[k]-cmin[k])*a;
}
}
}
//TagRealArray force = m->CreateTag("FORCE",DATA_REAL,CELL,NONE,3);
TagRealArray bc = m->CreateTag("BOUNDARY_CONDITION_VELOCITY",DATA_REAL,FACE,FACE,7);
......
#include "inmost.h"
using namespace INMOST;
//TODO:
//attach scale to rescale grid
//attach slice grid to cut cylinder
//option to set shear rate / pressure drop / maximum velocity
//setup BC
//cmake compilation
//setup on mesh:
//
// FORCE - 3 entries
//
// BOUNDARY_CONDITION_VELOCITY - 7 entries
// r = (a5,a6,a7)
// n.(a1*u + a2*t) = n.r
// (I - nn.)(a3*u + a4*t) = (I-nn.)r
//
// BOUNDARY_CONDITION_PRESSURE - 1 entry
//
// REFERENCE_VELOCITY - 3 entries
//
// REFERENCE_PRESSURE - 1 entry
//
double Len = 8;
double Diam = 1;
double shear = 10;//0; //shear rate
double dp = 0;//36e-6;
double vmax = 0;
typedef Storage::real real;
int main(int argc, char ** argv)
{
if( argc < 2 )
{
std::cout << "Usage: " << argv[0] << " mesh is2D=1 height=0.6 length=0.2 mesh_out=grid_out.pmf" << std::endl;
return 0;
}
Mesh * m = new Mesh;
m->SetFileOption("VERBOSITY","2");
try
{
m->Load(argv[1]);
}
catch(...)
{
std::cout << "Cannot load the mesh " << argv[1] << std::endl;
return -1;
}
std::string fout = "grid_out.pmf";
std::string
double cmax[3] = {-1.0e20,-1.0e20,-1.0e20}, cmin[3] = {1.0e20,1.0e20,1.0e20};
for(Mesh::iteratorNode n = m->BeginNode(); n != m->EndNode(); ++n)
{
Storage::real_array c = n->Coords();
for(int k = 0; k < 3; ++k)
{
if( cmax[k] < c[k] ) cmax[k] = c[k];
if( cmin[k] > c[k] ) cmin[k] = c[k];
}
}
Len = cmax[2] - cmin[2];
Diam = sqrt(pow(cmax[0] - cmin[0],2) + pow(cmax[1]-cmin[1],2))/sqrt(2);
if( argc > 5 ) fout = std::string(argv[5]);
//TagRealArray force = m->CreateTag("FORCE",DATA_REAL,CELL,NONE,3);
TagRealArray bc = m->CreateTag("BOUNDARY_CONDITION_VELOCITY",DATA_REAL,FACE,FACE,7);
TagReal bcp = m->CreateTag("BOUNDARY_CONDITION_PRESSURE",DATA_REAL,FACE,FACE,1);
TagRealArray uvw = m->CreateTag("UVW",DATA_REAL,CELL,NONE,3);
TagReal p = m->CreateTag("P",DATA_REAL,CELL,NONE,1);
//cut part of the grid
rMatrix n(3,1), x(3,1);
for(Mesh::iteratorCell it = m->BeginCell(); it != m->EndCell(); ++it)
{
it->Centroid(x.data());
if( x(0,0) < length && x(1,0) < height )
it->Delete();
}
//this should not be needed?
for(Mesh::iteratorFace it = m->BeginFace(); it != m->EndFace(); ++it) if( it->Boundary() )
it->FixNormalOrientation();
for(Mesh::iteratorCell it = m->BeginCell(); it != m->EndCell(); ++it)
{
uvw(*it,3,1).Zero();
p[*it] = 0;
}
if( is2D ) std::cout << "2D setup" << std::endl;
else std::cout << "3D setup" << std::endl;
std::cout << "height=" << height << " length=" << length << std::endl;
for(Mesh::iteratorFace it = m->BeginFace(); it != m->EndFace(); ++it) if( it->Boundary() )
{
it->Centroid(x.data());
it->UnitNormal(n.data());
if( is2D && (fabs(n(2,0)-1) < 1.0-3 || fabs(n(2,0)+1) < 1.0e-3 ) ) //roller bc
{
bc[*it][0] = 1;
bc[*it][1] = 0;
bc[*it][2] = 0;
bc[*it][3] = 1;
bc[*it][4] = 0;
bc[*it][5] = 0;
bc[*it][6] = 0;
}
else if( fabs(n(0,0)+1) < 1.0e-3 && x(0,0) < 1.0e-3 ) //inlet
{
bc[*it][0] = 1;
bc[*it][1] = 0;
bc[*it][2] = 1;
bc[*it][3] = 0;
bc[*it][4] = 1;
bc[*it][5] = 0;
bc[*it][6] = 0;
}
else if( fabs(n(0,0)-1) < 1.0e-3 ) //outlet
bcp[*it] = 0;
else //no-slip
{
bc[*it][0] = 1;
bc[*it][1] = 0;
bc[*it][2] = 1;
bc[*it][3] = 0;
bc[*it][4] = 0;
bc[*it][5] = 0;
bc[*it][6] = 0;
}
}
std::cout << "Saving output to " << fout << std::endl;
m->Save(fout);
return 0;
}
......@@ -12,9 +12,9 @@ namespace INMOST {
class SolverParameters;
enum SolverVerbosityLevel {
Level0 = 0,
Level1 = 1,
Level2 = 2
SolverVerbosityLevel0 = 0,
SolverVerbosityLevel1 = 1,
SolverVerbosityLevel2 = 2
};
/// Main class to set and solve linear system.
......
......@@ -17,7 +17,7 @@ namespace INMOST {
void Solver::SetParameterReal(std::string name, INMOST_DATA_REAL_TYPE value) { SetParameter(name, to_string(value)); }
Solver::Solver(std::string solverName, std::string prefix, INMOST_MPI_Comm _comm) :
versbosity(SolverVerbosityLevel::Level0), preconditioner_time(0), iterations_time(0), is_solved(false) {
versbosity(SolverVerbosityLevel0), preconditioner_time(0), iterations_time(0), is_solved(false) {
std::string lowerName = string_to_lower(solverName);
this->solver = SolverMaster::getSolver(solverName);
this->prefix = string_to_lower(prefix);
......@@ -164,7 +164,7 @@ namespace INMOST {
//INMOST::MPIBarrier();
iterations_time = Timer() - iterations_timer_start;
if (versbosity > SolverVerbosityLevel::Level1) {
if (versbosity > SolverVerbosityLevel1) {
#if defined(USE_MPI)
int rank;
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment