diff --git a/ConSol/ConSol/Mesh/BC/BCFAC.cpp b/ConSol/ConSol/Mesh/BC/BCFAC.cpp index c8115f5..8609a7e 100644 --- a/ConSol/ConSol/Mesh/BC/BCFAC.cpp +++ b/ConSol/ConSol/Mesh/BC/BCFAC.cpp @@ -1,30 +1,7 @@ // // Created by Fabian Moenkeberg on 2019-11-27. // #include "BCFAC.hpp" -void BCFAC::initialize(MeshBase *mesh, Configuration config) { - gamma = 1.4; - U_in.resize(4); - U_out.resize(4); - int M0 = 1; - double r_in = 1; - double u1_in = std::sqrt(gamma)*M0; - double u2_in =0; - double p_in = 1; - -// double p_out = 1.3; -// double r_out = r_in*(gamma-1+(gamma+1)*p_out)/(gamma+1+(gamma-1)*p_out); -// double u1_out = std::sqrt(gamma) + std::sqrt(2)*(1-p_out)/std::sqrt(gamma-1+p_out*(gamma+1)); -// double u2_out = 0; - - double E_in = p_in/(gamma - 1) + 0.5*(u1_in*u1_in+u2_in*u2_in)*r_in; -// double E_out = p_out/(gamma - 1) + 0.5*(u1_out*u1_out+u2_out*u2_out)*r_out; - - U_in = {r_in, r_in*u1_in, r_in*u2_in, E_in}; - U_out = {r_in, r_in*u1_in, r_in*u2_in, E_in}; -// U_in = U_out; -// U_out = U_in; -} \ No newline at end of file diff --git a/ConSol/ConSol/Mesh/BC/BCFAC.hpp b/ConSol/ConSol/Mesh/BC/BCFAC.hpp index 6c35218..fcfb912 100644 --- a/ConSol/ConSol/Mesh/BC/BCFAC.hpp +++ b/ConSol/ConSol/Mesh/BC/BCFAC.hpp @@ -1,21 +1,21 @@ // // Created by Fabian Moenkeberg on 2019-11-27. // #ifndef CONSOL_BCFAC_HPP #define CONSOL_BCFAC_HPP -#include "SlipWallBC.hpp" +#include "OpenFlowBC.hpp" -class BCFAC :public SlipWallBC{ +class BCFAC :public OpenFlowBC{ public: BCFAC(){ this->name = FAC; + alpha = 0; + M = 1; }; - - virtual void initialize(MeshBase *mesh, Configuration config); }; #endif //CONSOL_BCFAC_HPP diff --git a/ConSol/ConSol/Mesh/BC/SlipWallBC.cpp b/ConSol/ConSol/Mesh/BC/SlipWallBC.cpp index c52e754..85dffe8 100644 --- a/ConSol/ConSol/Mesh/BC/SlipWallBC.cpp +++ b/ConSol/ConSol/Mesh/BC/SlipWallBC.cpp @@ -1,220 +1,201 @@ // // Created by Fabian Moenkeberg on 2019-05-01. // #include "SlipWallBC.hpp" //#include "flux_farfield.cpp" void SlipWallBC::initialize(MeshBase *mesh, Configuration config) { gamma = 1.4; U_in.resize(4); U_out.resize(4); - double r_in = 1; - double u1_in = std::sqrt(gamma); - double u2_in =0; - double p_in = 1; - - double p_out = 1.3; - double r_out = r_in*(gamma-1+(gamma+1)*p_out)/(gamma+1+(gamma-1)*p_out); - double u1_out = std::sqrt(gamma) + std::sqrt(2)*(1-p_out)/std::sqrt(gamma-1+p_out*(gamma+1)); - double u2_out = 0; double E_in = p_in/(gamma - 1) + 0.5*(u1_in*u1_in+u2_in*u2_in)*r_in; double E_out = p_out/(gamma - 1) + 0.5*(u1_out*u1_out+u2_out*u2_out)*r_out; U_in = {r_in, r_in*u1_in, r_in*u2_in, E_in}; U_out = {r_out, r_out*u1_out, r_out*u2_out, E_out}; -// U_in = U_out; -// U_out = U_in; } void SlipWallBC::updateBoundary(std::vector > &U, MeshBase *mesh, double t){ -// std::cout << "testB0\n"; -// for (int k = 0; k < sysDim; k++){ -// std::cout << "testB"<< k<<"\n"; -// BC[k].updateBoundary(U, mesh); -// } -// int nCellsTot = mesh->getNcellsTot(); -//// std::cout << "testB1\n"; -// for (int i = 0; i < nCellsTot; i++){ -//// std::cout << "testBi\n"; -// U[i] = U[mesh->getBCMap(i)]; -// } - std::vector n; int nBC = mesh->getNbc(); for (int i = 0; i < nBC; i++){ U[mesh->getBCMap(3*i+0)] = U[mesh->getBCMap(3*i+1)]; int iEdge = mesh->getBCMap(3*i+2); if (iEdge == -1){ // internal hybrid update U[mesh->getBCMap(3 * i + 0)] = U[mesh->getBCMap(3 * i + 1)]; }else { std::vector edgesAtBC = mesh->getEdgesAtBC(iEdge); switch (edgesAtBC[1]) { case 2: { // slipwall BC -// std::vector edgesAtBC = mesh->getEdgesAtBC(mesh->getBCMap(3*i+2)); n = mesh->getNormal(edgesAtBC[0]); if (edgesAtBC[2] == 1) { //means left side of edge is out & normal is INWARDS pointing n[0] = -n[0]; n[1] = -n[1]; } flux_slipWall(U[mesh->getBCMap(3 * i + 1)], n, U[mesh->getBCMap(3 * i + 0)]); break; } case 0: { // means inflow BC -// std::vector edgesAtBC = mesh->getEdgesAtBC(mesh->getBCMap(3*i+2)); n = mesh->getNormal(edgesAtBC[0]); if (edgesAtBC[2] == 1) { //means left side of edge is out & normal is INWARDS pointing n[0] = -n[0]; n[1] = -n[1]; } U[mesh->getBCMap(3 * i + 0)] = flux_farfield(U[mesh->getBCMap(3 * i + 1)], n, U_in); break; } case 1: { // means outflow BC -// std::vector edgesAtBC = mesh->getEdgesAtBC(mesh->getBCMap(3*i+2)); n = mesh->getNormal(edgesAtBC[0]); if (edgesAtBC[2] == 1) { //means left side of edge is out & normal is INWARDS pointing n[0] = -n[0]; n[1] = -n[1]; } - +// U[mesh->getBCMap(3 * i + 0)] = flux_farfield(U[mesh->getBCMap(3 * i + 1)], n, U[mesh->getBCMap(3 * i + 1)]); U[mesh->getBCMap(3 * i + 0)] = flux_farfield(U[mesh->getBCMap(3 * i + 1)], n, U_out); +// U[mesh->getBCMap(3 * i + 0)] = U[mesh->getBCMap(3 * i + 1)]; break; } } } } int nBC2 = mesh->getNbc2()/2; int dimSys = mesh->getNelem(); std::vector tmp; std::vector tmp2; for (int i = 0; i < nBC2; i++){ tmp = U[mesh->getBCMap2(3*i + 1)]; tmp2 = U[mesh->getBCMap2(3*(i + nBC2) + 1)]; for (int j = 0; j < dimSys; j++){ U[mesh->getBCMap2(3*i + 0)][j] = 0.5*(tmp[j] + tmp2[j]); } } } void SlipWallBC::updateEdges(std::vector>>> &UReconstr, MeshBase *mesh, double t){ std::vector> edgesAtBC = mesh->getEdgesAtBC(); int nelem = mesh->getNelem(); int nQuad = UReconstr[0][0].size(); int nEdgesBC = edgesAtBC.size(); std::vector n; for(int i = 0; i < nEdgesBC; i++) { nQuad = UReconstr[0][edgesAtBC[i][0]].size(); for (int k = 0; k < nQuad; k++){ -// if ((edgesAtBC[i][1] == 1)||(edgesAtBC[i][1] == 3)){// top and bottom wall with slipwall BC if ((edgesAtBC[i][1] == 2)){// slipwall BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); flux_slipWall(UReconstr[1][edgesAtBC[i][0]][k], n, UReconstr[0][edgesAtBC[i][0]][k]); }else{//means left side of edge is out & normal is INWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); n[0] = -n[0]; n[1] = -n[1]; flux_slipWall(UReconstr[0][edgesAtBC[i][0]][k], n, UReconstr[1][edgesAtBC[i][0]][k]); } }else if ((edgesAtBC[i][1] == 0)){ // means inflow BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); UReconstr[0][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[1][edgesAtBC[i][0]][k], n, U_in); }else{//means left side of edge is out & normal is INWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); n[0] = -n[0]; n[1] = -n[1]; UReconstr[1][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[0][edgesAtBC[i][0]][k], n, U_in); } -// }else{ //right wall, means outflow BC }else if ((edgesAtBC[i][1] == 1)){ // means outflow BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); UReconstr[0][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[1][edgesAtBC[i][0]][k], n, U_out); +// UReconstr[0][edgesAtBC[i][0]][k] = UReconstr[1][edgesAtBC[i][0]][k]; +// UReconstr[0][edgesAtBC[i][0]][k] = U_out; }else{//means left side of edge is out & normal is INWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); n[0] = -n[0]; n[1] = -n[1]; UReconstr[1][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[0][edgesAtBC[i][0]][k], n, U_out); +// UReconstr[1][edgesAtBC[i][0]][k] = UReconstr[0][edgesAtBC[i][0]][k]; +// UReconstr[1][edgesAtBC[i][0]][k] = U_out; } } } if ((edgesAtBC[i][1] == 11)){ // means hybrid Connection BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing for (int k = 0; k < nQuad; k++){ UReconstr[0][edgesAtBC[i][0]][k] = UReconstr[0][edgesAtBC[i][3]][0]; } }else{//means left side of edge is out & normal is INWARDS pointing for ( int l = 0; l < nelem; l++){ UReconstr[1][edgesAtBC[i][0]][0][l] = 0; for (int k = 0; k < nQuad; k++){ UReconstr[1][edgesAtBC[i][0]][0][l] += UReconstr[1][edgesAtBC[i][3]][k][l]/nQuad; } } } }else if ((edgesAtBC[i][1] == 14)){ // means hybrid QUADQUAD 2 QUAD or QUAD 2 QUADQUAD Connection BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing UReconstr[0][edgesAtBC[i][0]]= UReconstr[1][edgesAtBC[i][3]]; }else{//means left side of edge is out & normal is INWARDS pointing UReconstr[1][edgesAtBC[i][0]]= UReconstr[0][edgesAtBC[i][3]]; } } } } void SlipWallBC::updateEdges(std::vector>>> &UReconstr, MeshBase *mesh, double t, int i){ std::vector> edgesAtBC = mesh->getEdgesAtBC(); int nQuad = UReconstr[0][0].size(); int nEdgesBC = edgesAtBC.size(); std::vector n; for (int k = 0; k < nQuad; k++){ if ((edgesAtBC[i][1] == 2)){// slipwall BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); flux_slipWall(UReconstr[1][edgesAtBC[i][0]][k], n, UReconstr[0][edgesAtBC[i][0]][k]); }else{//means left side of edge is out & normal is INWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); n[0] = -n[0]; n[1] = -n[1]; flux_slipWall(UReconstr[0][edgesAtBC[i][0]][k], n, UReconstr[1][edgesAtBC[i][0]][k]); } }else if ((edgesAtBC[i][1] == 0)){ // means inflow BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); UReconstr[0][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[1][edgesAtBC[i][0]][k], n, U_in); }else{//means left side of edge is out & normal is INWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); n[0] = -n[0]; n[1] = -n[1]; UReconstr[1][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[0][edgesAtBC[i][0]][k], n, U_in); } }else if ((edgesAtBC[i][1] == 1)){ // means outflow BC if (edgesAtBC[i][2] == 0) { //means right side of edge is out & normal is OUTWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); - UReconstr[0][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[1][edgesAtBC[i][0]][k], n, U_out); +// UReconstr[0][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[1][edgesAtBC[i][0]][k], n, U_out); +// UReconstr[0][edgesAtBC[i][0]][k] = UReconstr[1][edgesAtBC[i][0]][k]; + UReconstr[0][edgesAtBC[i][0]][k] = U_out; }else{//means left side of edge is out & normal is INWARDS pointing n = mesh->getNormal(edgesAtBC[i][0]); n[0] = -n[0]; n[1] = -n[1]; - UReconstr[1][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[0][edgesAtBC[i][0]][k], n, U_out); +// UReconstr[1][edgesAtBC[i][0]][k] = flux_farfield(UReconstr[0][edgesAtBC[i][0]][k], n, U_out); +// UReconstr[1][edgesAtBC[i][0]][k] = UReconstr[0][edgesAtBC[i][0]][k]; + UReconstr[1][edgesAtBC[i][0]][k] = U_out; } } } } void SlipWallBC::flux_slipWall(std::vector U, std::vector n, std::vector &ret){ ret = U; std::vector V_in{U[1]/U[0],U[2]/U[0]}; double vn = V_in[0]*n[0] + V_in[1]*n[1]; ret[1] = U[0]*(V_in[0]-2*vn*n[0]); ret[2] = U[0]*(V_in[1]-2*vn*n[1]); } \ No newline at end of file diff --git a/ConSol/ConSol/Mesh/BC/SlipWallBC.hpp b/ConSol/ConSol/Mesh/BC/SlipWallBC.hpp index ee7171d..f049cc8 100644 --- a/ConSol/ConSol/Mesh/BC/SlipWallBC.hpp +++ b/ConSol/ConSol/Mesh/BC/SlipWallBC.hpp @@ -1,54 +1,66 @@ // // Created by Fabian Moenkeberg on 2019-05-01. // #ifndef CONSOL_SLIPWALLBC_HPP #define CONSOL_SLIPWALLBC_HPP #include #include "BcBase.hpp" #include "../../Base/Configuration.hpp" #include "../MeshBase.hpp" class SlipWallBC : public BcBase{ protected: int sysDim; double gamma = 1.4; + double r_in; + double u1_in; + double u2_in; + double p_in ; + + double p_out; + double r_out; + double u1_out; + double u2_out; + std::vector BC; std::vector U_in; std::vector U_out; public: SlipWallBC(){ this->name = SLIPWALL; + r_in = 1; + u1_in = std::sqrt(gamma); + u2_in =0; + p_in = 1; + + p_out = 1.3; + r_out = r_in*(gamma-1+(gamma+1)*p_out)/(gamma+1+(gamma-1)*p_out); + u1_out = std::sqrt(gamma) + std::sqrt(2)*(1-p_out)/std::sqrt(gamma-1+p_out*(gamma+1)); + u2_out = 0; }; -// MixedBC_2D(int sysdim){ -// this->name = MIXEDBC; -// BC.resize(sysdim); -// for (int i = 0; i < sysdim; i++){ -// BC[i] = MixedBcComponent_2D(); -// } -// }; SlipWallBC(std::vector BC){ this->sysDim = (int) BC.size(); this->BC = BC; this->name = SLIPWALL; } ~SlipWallBC(){} virtual void initialize(MeshBase *mesh, Configuration config); virtual void updateBoundary(std::vector> &U, MeshBase *mesh, double t); virtual void updateEdges(std::vector>>> &UReconstr, MeshBase *mesh, double t); virtual void updateEdges(std::vector>>> &UReconstr, MeshBase *mesh, double t, int i); void flux_slipWall(std::vector U, std::vector n, std::vector &ret); }; #endif //CONSOL_SLIPWALLBC_HPP diff --git a/ConSol/Examples/HybridTest/FlowAroundCylinder/FAC8QUAD.geo b/ConSol/Examples/HybridTest/FlowAroundCylinder/FAC8QUAD.geo new file mode 100644 index 0000000..c7f49c7 --- /dev/null +++ b/ConSol/Examples/HybridTest/FlowAroundCylinder/FAC8QUAD.geo @@ -0,0 +1,269 @@ +// Variables +boxdimX = 3; +boxdimY = 3; + + +N = 4*28*3; +M = N; +N0 = 2; +gridsizeX = boxdimX/(N); +gridsizeY = boxdimY/(M); + +If (gridsizeX < gridsizeY) + gridsize = gridsizeX; +Else + gridsize = gridsizeY; +EndIf + +If (boxdimX < boxdimY) + boxdim = boxdimX; +Else + boxdim = boxdimY; +EndIf + +Nx1 = (N)/3 - N0; +Nx2 = ((N)/3); +Nx3 = Nx1*2 + N0; +Ny1 = (M)/3 - N0; +Ny2 = ((M)/3); +Ny3 = Ny1; + +X0 = 0; +X1 = X0 + (Nx1+N0)*gridsizeX; +X2 = X1 + Nx2*gridsizeX; +X3 = X2 + (Nx3+N0)*gridsizeX; + +Y0 = 0; +Y1 = Y0 + (Ny1+N0)*gridsizeY; +Y2 = Y1 + Ny2*gridsizeY; +Y3 = Y2 + (Ny3+N0)*gridsizeY; + +D = (X2-X1)*7/8/2; + +/*Nx1 = (N-1)*(X1-X0)/(X3-X0)+1 - N0; +Nx2 = ((N-1)*(X2-X1)/(X3-X0)) + 1; +Nx3 = (N-1)*(X3-X2)/(X3-X0)+1 - N0; +Ny1 = (N-1)*(Y1-Y0)/(Y3-Y0)+1 - N0; +Ny2 = ((N-1)*(Y2-Y1)/(Y3-Y0)) + 1; +Ny3 = (N-1)*(Y3-Y2)/(Y3-Y0)+1 - N0;*/ + +// Points +Point(1) = {X0,Y0,0, gridsize}; +Point(2) = {X3,Y0,0, gridsize*3}; +Point(3) = {X3,Y3,0, gridsize*3}; +Point(4) = {X0,Y3,0, gridsize}; + +Point(5) = {X1-N0*gridsize,Y0,0, gridsize}; +Point(6) = {X1,Y0,0, gridsize}; +Point(7) = {X2,Y0,0, gridsize}; +Point(8) = {X2+N0*gridsize,Y0,0, gridsize}; + +Point(9) = {X3, Y1-N0*gridsize,0, gridsize}; +Point(10) = {X3, Y1,0, gridsize}; +Point(11) = {X3, Y2,0, gridsize}; +Point(12) = {X3, Y2+N0*gridsize,0, gridsize}; + +Point(13) = {X2+N0*gridsize,Y3,0, gridsize}; +Point(14) = {X2,Y3,0, gridsize}; +Point(15) = {X1,Y3,0, gridsize}; +Point(16) = {X1-N0*gridsize,Y3,0, gridsize}; + +Point(17) = {X0, Y2+N0*gridsize,0, gridsize}; +Point(18) = {X0, Y2,0, gridsize}; +Point(19) = {X0, Y1,0, gridsize}; +Point(20) = {X0, Y1-N0*gridsize,0, gridsize}; + +Point(21) = {X1-N0*gridsize,Y1-N0*gridsize,0, gridsize}; +Point(22) = {X1,Y1-N0*gridsize,0, gridsize}; +Point(23) = {X1,Y1,0, gridsize}; +Point(24) = {X1-N0*gridsize,Y1,0, gridsize}; + +Point(25) = {X2,Y1-N0*gridsize,0, gridsize}; +Point(26) = {X2+N0*gridsize,Y1-N0*gridsize,0, gridsize}; +Point(27) = {X2+N0*gridsize,Y1,0, gridsize}; +Point(28) = {X2,Y1,0, gridsize}; + +Point(29) = {X2,Y2,0, gridsize}; +Point(30) = {X2+N0*gridsize,Y2,0, gridsize}; +Point(31) = {X2+N0*gridsize,Y2+N0*gridsize,0, gridsize}; +Point(32) = {X2,Y2+N0*gridsize,0, gridsize}; + +Point(33) = {X1-N0*gridsize,Y2,0, gridsize}; +Point(34) = {X1,Y2,0, gridsize}; +Point(35) = {X1,Y2+N0*gridsize,0, gridsize}; +Point(36) = {X1-N0*gridsize,Y2+N0*gridsize,0, gridsize}; + +//Circle +Point(37) = {(X1+X2)/2,(Y0+Y3)/2 - D,0,1.0}; +Point(38) = {(X1+X2)/2 + D,(Y0+Y3)/2,0,1.0}; +Point(39) = {(X1+X2)/2,(Y0+Y3)/2+D,0,1.0}; +Point(40) = {(X1+X2)/2 - D,(Y0+Y3)/2,0,1.0}; +Point(41) = {(X1+X2)/2,(Y0+Y3)/2, 0 ,1.0}; + + +//Lines +Line(1) = {1,5}; +Line(2) = {5,6}; +Line(3) = {6,7}; +Line(4) = {7,8}; +Line(5) = {8,2}; +Line(6) = {2,9}; +Line(7) = {9,10}; +Line(8) = {10,11}; +Line(9) = {11,12}; +Line(10) = {12,3}; +Line(11) = {13,3}; +Line(12) = {14,13}; +Line(13) = {15,14}; +Line(14) = {16,15}; +Line(15) = {4,16}; +Line(16) = {17,4}; +Line(17) = {18,17}; +Line(18) = {19,18}; +Line(19) = {20,19}; +Line(20) = {1,20}; +Line(21) = {20,21}; +Line(22) = {21,22}; +Line(23) = {22,25}; +Line(24) = {25,26}; +Line(25) = {26,9}; +Line(26) = {19,24}; +Line(27) = {24,23}; +Line(28) = {23,28}; +Line(29) = {28,27}; +Line(30) = {27,10}; +Line(31) = {18,33}; +Line(32) = {33,34}; +Line(33) = {34,29}; +Line(34) = {29,30}; +Line(35) = {30,11}; +Line(36) = {17,36}; +Line(37) = {36,35}; +Line(38) = {35,32}; +Line(39) = {32,31}; +Line(40) = {31,12}; +Line(41) = {5,21}; +Line(42) = {21,24}; +Line(43) = {24,33}; +Line(44) = {33,36}; +Line(45) = {36,16}; +Line(46) = {6,22}; +Line(47) = {22,23}; +Line(48) = {23,34}; +Line(49) = {34,35}; +Line(50) = {35,15}; +Line(51) = {7,25}; +Line(52) = {25,28}; +Line(53) = {28,29}; +Line(54) = {29,32}; +Line(55) = {32,14}; +Line(56) = {8,26}; +Line(57) = {26,27}; +Line(58) = {27,30}; +Line(59) = {30,31}; +Line(60) = {31,13}; +Circle(61) = {37, 41, 38}; +Circle(62) = {38, 41, 39}; +Circle(63) = {39, 41, 40}; +Circle(64) = {40, 41, 37}; + +Line Loop(1) = {1,41,-21,-20}; +Line Loop(2) = {2,46,-22,-41}; +Line Loop(3) = {3,51,-23,-46}; +Line Loop(4) = {4,56,-24,-51}; +Line Loop(5) = {5,6,-25,-56}; +Line Loop(6) = {21,42,-26,-19}; +Line Loop(7) = {22,47,-27,-42}; +Line Loop(8) = {23,52,-28,-47}; +Line Loop(9) = {24,57,-29,-52}; +Line Loop(10) = {25,7,-30,-57}; +Line Loop(11) = {26,43,-31,-18}; +Line Loop(12) = {27,48,-32,-43}; +Line Loop(13) = {28,53,-33,-48}; +Line Loop(14) = {29,58,-34,-53}; +Line Loop(15) = {30,8,-35,-58}; +Line Loop(16) = {31,44,-36,-17}; +Line Loop(17) = {32,49,-37,-44}; +Line Loop(18) = {33,54,-38,-49}; +Line Loop(19) = {34,59,-39,-54}; +Line Loop(20) = {35,9,-40,-59}; +Line Loop(21) = {36,45,-15,-16}; +Line Loop(22) = {37,50,-14,-45}; +Line Loop(23) = {38,55,-13,-50}; +Line Loop(24) = {39,60,-12,-55}; +Line Loop(25) = {40,10,-11,-60}; + +Line Loop(26) = {61,62,63,64}; + +Plane Surface(1) = {1}; +Plane Surface(2) = {2}; +Plane Surface(3) = {3}; +Plane Surface(4) = {4}; +Plane Surface(5) = {5}; +Plane Surface(6) = {6}; +Plane Surface(7) = {7}; +Plane Surface(8) = {8}; +Plane Surface(9) = {9}; +Plane Surface(10) = {10}; +Plane Surface(11) = {11}; +Plane Surface(12) = {12}; +Plane Surface(13) = {13,-26}; +Plane Surface(14) = {14}; +Plane Surface(15) = {15}; +Plane Surface(16) = {16}; +Plane Surface(17) = {17}; +Plane Surface(18) = {18}; +Plane Surface(19) = {19}; +Plane Surface(20) = {20}; +Plane Surface(21) = {21}; +Plane Surface(22) = {22}; +Plane Surface(23) = {23}; +Plane Surface(24) = {24}; +Plane Surface(25) = {25}; + +Physical Surface("tri_0") = {13}; +Physical Surface("qdr_0") = {1}; +Physical Surface("qdr_1") = {3}; +Physical Surface("qdr_2") = {5}; +Physical Surface("qdr_3") = {11}; +Physical Surface("qdr_4") = {15}; +Physical Surface("qdr_5") = {21}; +Physical Surface("qdr_6") = {23}; +Physical Surface("qdr_7") = {25}; +Physical Surface("mid_0") = {8}; +Physical Surface("mid_1") = {12}; +Physical Surface("mid_2") = {14}; +Physical Surface("mid_3") = {18}; +Physical Surface("q2q_0") = {2}; +Physical Surface("q2q_1") = {4}; +Physical Surface("q2q_2") = {6}; +Physical Surface("q2q_3") = {10}; +Physical Surface("q2q_4") = {16}; +Physical Surface("q2q_5") = {20}; +Physical Surface("q2q_6") = {22}; +Physical Surface("q2q_7") = {24}; +Physical Surface("rst_1") = {7}; +Physical Surface("rst_2") = {9}; +Physical Surface("rst_3") = {17}; +Physical Surface("rst_4") = {19}; + +Physical Line("inflow") = {20,19,18,17,16}; +Physical Line("outflow") = {6,7,8,9,10}; +Physical Line("slipWall") = {1,2,3,4,5,15,14,13,12,11,61,62,63,64}; +//Physical Line("slipWall") = {1,2,3,4,5,15,14,13,12,11}; +Physical Line("q2i") = {23,58,38,43,22,47,27,42,24,57,29,52,34,59,39,54,32,49,37,44}; +Physical Line("t2i") = {28,53,33,48}; +Physical Line("q2q") = {41,46,51,56,45,50,55,60,21,26,31,36,25,30,35,40}; + +Transfinite Curve{1,21,26,31,36,15} = Nx1+1; +Transfinite Curve{5,25,30,35,40,11} = Nx3+1; +Transfinite Curve{20,41,46,51,56,6} = Ny1+1; +Transfinite Curve{16,45,50,55,60,10} = Ny3+1; +Transfinite Curve{2,22,27,32,37,14,4,24,29,34,39,12,19,42,47,52,57,7,17,44,49,54,59,9} = N0+1; +Transfinite Curve{3,23,28,33,38,13} = Nx2+1; +Transfinite Curve{18,43,48,53,58,8} = Ny2+1; + +Transfinite Curve {61,62,63,64} = 2*D*N/boxdim/1.2+1 Using Progression 1; + +Transfinite Surface{1,2,3,4,5,6,7,8,9,10,11,12,14,15,16,17,18,19,20,21,22,23,24,25}; +Recombine Surface{1,2,3,4,5,6,7,8,9,10,11,12,14,15,16,17,18,19,20,21,22,23,24,25}; diff --git a/ConSol/Examples/HybridTest/FlowAroundCylinder/paramFAC.in b/ConSol/Examples/HybridTest/FlowAroundCylinder/paramFAC.in index 904699c..a859492 100644 --- a/ConSol/Examples/HybridTest/FlowAroundCylinder/paramFAC.in +++ b/ConSol/Examples/HybridTest/FlowAroundCylinder/paramFAC.in @@ -1,106 +1,106 @@ grid { load_grid true dim 2 x_min 0 x_max 1 Nx_cell { 10 } nghostcells 2 - grid_path ../../Examples/HybridTest/FlowAroundCylinder/FAC8QUAD_0.msh + grid_path ../../Examples/HybridTest/FlowAroundCylinder/FAC8QUAD.msh type hybrid } numeric { - nthreads 4 + nthreads 0 time_scheme FE - cfl 0.5 + cfl 0.8 max_iter 50000 final_time 0.35 min_residue 1.0e-6 reconstruct weno stencilChoice otf order 3 degPoly 2 interpolationBase Multiquadratics shapeParameter 0.1 reconstrOnCharact false evaluationMethod meanvalue fluxScheme LaxFried limiter noLimiter tInclude { } } model { model Euler source noSource } constants { dl 0.445 ul 0.698 vl 0.0 pl 3.528 dr 0.5 ur 0.0 vr 0.0 pr 0.571} initial_condition { y initFAC } exact_soln { available no y } boundary_condition { bc FAC xleft { type slipWall } xright { type slipWall } yleft { type slipWall } yright { type slipWall } } output { output_path ../../Examples/HybridTest/FlowAroundCylinder } import { import false import_path ../../Examples/HybridTest/FlowAroundCylinder/solutionref1.h5 } material { gamma 1.4 gas_const 1.0 viscosity { model constant mu_ref 0.0 } prandtl 0.72 }