diff --git a/hemoFlow.cpp b/hemoFlow.cpp index 6cd5193..e8b9470 100644 --- a/hemoFlow.cpp +++ b/hemoFlow.cpp @@ -48,6 +48,7 @@ string workingFolder; T simLength; T saveFreqTime; T checkpointFreqTime; +T convergenceTolerance = 1e-11; // Vector of openings vector openings; @@ -75,7 +76,7 @@ int findIndex(const unsigned short *array, int arraySize, unsigned short itemToF // *** Calculating LB parameters using Re on the inlet: Re = U_avg * D / nu void calcSimulationParameters(SimPar &sim, T dx, T dt = -1, T U_max_LB_ = 0.1) -{ +{ sim.C_l = dx; sim.U_max_lb = U_max_LB_; @@ -85,7 +86,7 @@ void calcSimulationParameters(SimPar &sim, T dx, T dt = -1, T U_max_LB_ = 0.1) nuInf_lb = nuInf * sim.C_t / sim.C_l / sim.C_l; T tau = 3.0*nuInf_lb+0.5; sim.omega = 1.0 / tau; - } + } else { sim.omega = 1.0; nuInf_lb = 0.5 / 3.0; @@ -187,7 +188,7 @@ int main(int argc, char *argv[]) << (std::string)global::argv(0) << " parameter-input-file.xml [-r]" << std::endl; return -1; } - + // Check if we have the checkpoint flag string checkpointFlag; bool isCheckpointed = false; @@ -216,17 +217,17 @@ int main(int argc, char *argv[]) workingFolder = "."; else workingFolder = paramXmlFileName.substr(0, folderIdx); - + // *** Load in data files try { pcout << "Loading in data file..." << std::endl; xml["simulation"]["outputDir"].read(outputFolder); - + outDir = workingFolder + "/" + outputFolder; // Check if output dir exists and accessible on the main processor - if(global::mpi().isMainProcessor()) + if(global::mpi().isMainProcessor()) if (dirExists(outDir) <= 0) { pcout << "Output folder " << outDir << " does not exist! Creating it...." << std::endl; mkpath(outDir.c_str(), 0777); @@ -243,7 +244,7 @@ int main(int argc, char *argv[]) xml["simulation"]["blockSize"].read(blockSize); xml["simulation"]["simLength"].read(simLength); xml["simulation"]["saveFrequency"].read(saveFreqTime); - + // Check for optional checkpoint argument try { xml["simulation"]["checkpointFrequency"].read(checkpointFreqTime); @@ -253,13 +254,22 @@ int main(int argc, char *argv[]) useCheckpoint = false; } + // Reading optional convergence tolerance + try { + xml["simulation"]["convergenceTolerance"].read(convergenceTolerance); + } + catch (PlbIOException& exception) { + pcout << "Warning: convergenceTolerance tag was not found in config. Using default value." << std::endl; + pcout << "Default convergence tolerance: " << convergenceTolerance << std::endl; + } + // Loading the input file string npzFileName; xml["geometry"]["file"].read(npzFileName); cnpy::npz_t geom_npz = cnpy::npz_load(workingFolder + "/" + npzFileName); pcout << "Input data elements: " << geom_npz.size() << std::endl; - - for (auto const& array : geom_npz) + + for (auto const& array : geom_npz) pcout << "->" << array.first << std::endl; // TODO - Consider hdf5-based input files (see the Turbulence branch) @@ -279,7 +289,7 @@ int main(int argc, char *argv[]) // Loading the geometry of the stent (if there is one) stentFlag = geom_npz["stent"]; - + if(stentFlag.shape.size() > 1) { // Check if there is data on the FD pcout << "Found flow diverter information to load." << std::endl; sfData = stentFlag.data(); @@ -297,22 +307,22 @@ int main(int argc, char *argv[]) // Calcualte simulation parameters here! pcout << "Setting LBM parameters..." << std::endl; - calcSimulationParameters(sim, sim_dx, sim_dt); + calcSimulationParameters(sim, sim_dx, sim_dt); // **** Processing openings **** pcout << "Processing openings..." << std::endl; - - // Loading information on openings + + // Loading information on openings cnpy::NpyArray openingIndex = geom_npz["openingIndex"]; auto* oiData = openingIndex.data(); - + cnpy::NpyArray openingRadius = geom_npz["openingRadius"]; auto* orData = openingRadius.data(); - + /* -- Not needed atm. cnpy::NpyArray openingQRatio = geom_npz["openingNormalizedQRatio"]; double* oqData = openingQRatio.data(); - + cnpy::NpyArray openingCenter = geom_npz["openingCenter"]; double* ocData = openingCenter.data(); */ @@ -326,17 +336,17 @@ int main(int argc, char *argv[]) setToFunction(*flagMatrix, flagMatrix->getBoundingBox(), FlagMaskDomain3D(gfData, 1)); pcout << "Creating sparse representation ..." << endl; - + //Create sparse representation MultiBlockManagement3D sparseBlockManagement = computeSparseManagement(*plb::reparallelize(*flagMatrix, blockSize, blockSize, blockSize), envelopeWidth); - + #if LES lattice = new MultiBlockLattice3D (sparseBlockManagement, defaultMultiBlockPolicy3D().getBlockCommunicator(), defaultMultiBlockPolicy3D().getCombinedStatistics(), defaultMultiBlockPolicy3D().getMultiCellAccess(), - new SmagorinskyRegularizedDynamics(sim.omega, cSmago)); + new SmagorinskyRegularizedDynamics(sim.omega, cSmago)); #else lattice = new MultiBlockLattice3D (sparseBlockManagement, defaultMultiBlockPolicy3D().getBlockCommunicator(), @@ -359,7 +369,7 @@ int main(int argc, char *argv[]) pcout << getMultiBlockInfo(*lattice) << endl; - // Loop through the openings in the datafile + // Loop through the openings in the datafile unsigned int numOpenings = openingRadius.shape[0]; pcout << "Number of openings in geometry: " << numOpenings << std::endl; pcout << "Opening flags: "; @@ -371,28 +381,28 @@ int main(int argc, char *argv[]) // Loop through the openings following the xml config for(unsigned int o=0; o < numOpenings; o++){ pcout << "Processing opening: " << o << std::endl; - + string xmlTagOpening = "opening_"+std::to_string(o); - + string name; xml["geometry"][xmlTagOpening]["name"].read(name); int type; xml["geometry"][xmlTagOpening]["type"].read(type); int label; xml["geometry"][xmlTagOpening]["label"].read(label); - int openingIdx = findIndex(oiData, numOpenings, label); - if(openingIdx == -1) + int openingIdx = findIndex(oiData, numOpenings, label); + if(openingIdx == -1) pcout << "ERROR: Opening label " << label << " was found in the config xml, but not in the geometry file!" << endl; - // Get the direction of the opening + // Get the direction of the opening int s = openingNormal.shape[1]; vec3d dir(onData[gT2D(s,openingIdx,0)], onData[gT2D(s,openingIdx,1)], onData[gT2D(s,openingIdx,2)]); // Create the opening auto *opening = new OpeningHandler(gfData, static_cast(label), static_cast(type), orData[openingIdx] / sim.C_l, dir); - + opening->setName(name); opening->setBCType(lattice); - string parameterStr; + string parameterStr; double parameter; xml["geometry"][xmlTagOpening]["parameter"].read(parameterStr); if(!parameterStr.empty()) @@ -410,20 +420,20 @@ int main(int argc, char *argv[]) if(type == OPENING_VELOCITY || type == OPENING_MURRAY || type == OUTLET_FREEFLOW){ opening->createPoiseauilleProfile(); // Normalized to max_vel = 1.0 (LBM units) // opening->normalizeFlowRate(); // Normalize to Q=1 (LBM units) - } + } else if(type == OPENING_PRESSURE) { opening->createConstantPressureProfile(); // Contant pressure = 1.0 (LBM density) } opening->printOpeningDetails(sim); - + openings.push_back(opening); } // Sanity check if(numOpenings != openings.size()) { - pcout << "**WARNING** The number of opening definitions don't match between the config and the geometry file! So how many openings do we actually have?" << endl; + pcout << "**WARNING** The number of opening definitions don't match between the config and the geometry file! So how many openings do we actually have?" << endl; } } @@ -432,7 +442,7 @@ int main(int argc, char *argv[]) << ": " << exception.what() << std::endl; return -1; } - + pcout << "*********** Simulation parameters *********** " << endl << "size [LU]: " << Nx << "x" << Ny << "x" << Nz << endl << "dx [m]: " << sim.C_l << endl @@ -440,7 +450,7 @@ int main(int argc, char *argv[]) << "dm [kg]: " << sim.C_m << endl << "omega: " << sim.omega << endl << "nu: " << 1./3. * (1./sim.omega - 0.5) << endl - << "U_max [lbm]: " << sim.U_max_lb << endl + << "U_max [lbm]: " << sim.U_max_lb << endl << "U_max [m/s]: " << sim.U_max_lb * sim.C_l / sim.C_t << endl << endl; int saveFrequency; @@ -448,17 +458,17 @@ int main(int argc, char *argv[]) pcout << "Saving frequency set to every " << saveFreqTime << " s (" << saveFrequency << " steps)." << endl; int checkpointFrequency = (int)round(checkpointFreqTime/sim.C_t); - - if(useCheckpoint) + + if(useCheckpoint) pcout << "Chekpointing will happen every " << checkpointFreqTime << " s (" << checkpointFrequency << " steps)." << endl; - + // Checkpoint file names relative to the output folder string chkParamFile = outDir+"/checkpoint_parameters.dat"; string chkDataFile = outDir+"/checkpoint_lattice.dat"; string chkParamFileOld = outDir+"/checkpoint_parameters_old.dat"; string chkDataFileOld = outDir+"/checkpoint_lattice_old.dat"; - // If there is data on porosity, set up porous layer in the simulation + // If there is data on porosity, set up porous layer in the simulation if(sfData != nullptr) { pcout << "Setting up porous layer for flow diverter..." << std::endl; porosityField = defaultGenerateMultiNTensorField3D(lattice->getMultiBlockManagement(), 1).release(); @@ -471,7 +481,7 @@ int main(int argc, char *argv[]) defineDynamics(*lattice, lattice->getBoundingBox(), new FlagMaskSingleDomain3D(gfData, 1), new BounceBack(1.0)); // TODO: add some reparallelize here, check if it plays nice with checkpointing - + pcout << "Initializing lattice in equilibrium..." << std::endl; initializeAtEquilibrium (*lattice, lattice->getBoundingBox(), 1.0, Array((T)0.,(T)0.,(T)0.) ); @@ -483,11 +493,11 @@ int main(int argc, char *argv[]) // iteration counter int stat_cycle = 0; - + // Check if the simulation was checkpointed if(isCheckpointed) { pcout << endl << "*********** Restoring checkpoint ***********" << endl; - + // Load in the iteration counter string chkParamFile = outDir+"/checkpoint_parameters.dat"; plb_ifstream ifile(chkParamFile.c_str()); @@ -499,18 +509,18 @@ int main(int argc, char *argv[]) pcout << "ERROR reading from the checkpoint parameter file: checkpoint_parameters.dat! Exiting..." << std::endl; return -1; } - + // Load the lattice loadBinaryBlock(*lattice, outDir+"/checkpoint_lattice.dat"); - + pcout << "Checkpoint at iteration " << stat_cycle << " loaded succesfully." << std::endl; } else { // If not, then let's chek the initial state and do a warm up. pcout << endl << "*********** Entering stationary warmup phase ***********" << endl; - + int convergenceSteps = 10*max(max(Nx, Ny), Nz); - T minDE = 1e-11; T dE = 100; T prevE = 0; - + T minDE = convergenceTolerance; T dE = 100; T prevE = 0; + // imposeOpenings(0.0); T cE = computeAverageEnergy(*lattice); @@ -529,38 +539,38 @@ int main(int argc, char *argv[]) while(abs(dE) > minDE && stat_cycle < convergenceSteps ) { lattice->collideAndStream(); - + T cE = computeAverageEnergy(*lattice); dE = cE - prevE; prevE = cE; - + if(stat_cycle % 500 == 0) { pcout << "Delta energy: " << abs(dE) << "/" << minDE << " Cycle: [" << stat_cycle << "/" << convergenceSteps <<"]" << std::endl; } - + stat_cycle++; } pcout << "Delta energy: " << abs(dE) << "/" << minDE << " Cycle: [" << stat_cycle << "/" << convergenceSteps <<"]" << std::endl; - + pcout << endl << "*********** Entering transient simulation phase ***********" << endl; - + pcout << "Saving time step 0..." << endl; // writeVTK(*lattice, 0); // writeNPZ(*lattice, 0); writeHDF5(*lattice, sim, 0, outDir); - + // Set the counter back stat_cycle = 0; } - + pcout << "Starting computation..." << endl; while(stat_cycle*sim.C_t <= simLength + sim.C_t) { - + if(stat_cycle % 200 == 0) { T cE = computeAverageEnergy(*lattice); pcout << "\rTime: " << stat_cycle*sim.C_t << "s / " << simLength << "s" << " [" << stat_cycle << " / " << std::round(simLength/sim.C_t) << "] " << " - Energy: " << cE << endl; - + // Capture numerical divergence if appears if (std::isnan(cE)){ pcout << "ERROR: NaN average energy! Saving state and stopping simulation" << std::endl; @@ -582,31 +592,31 @@ int main(int argc, char *argv[]) // Advance time stat_cycle++; - // Save output + // Save output if(stat_cycle % saveFrequency == 0) { pcout << "Writing output at: " << stat_cycle << " (" << stat_cycle*sim.C_t << " s)." << endl; // writeVTK(*lattice, sim, stat_cycle); // writeNPZ(*lattice, sim, stat_cycle); - writeHDF5(*lattice, sim, stat_cycle, outDir, porosityField); + writeHDF5(*lattice, sim, stat_cycle, outDir, porosityField); } - + if(useCheckpoint && (stat_cycle % checkpointFrequency == 0)) { // Overwriting previous checkpoint. Note: if failure happens during saving the checkpoint we cannot recover: TODO two step checkpoint if(global::mpi().isMainProcessor()) { if(fileExists(chkDataFile)){ - // Remove prev-previous checkpoint + // Remove prev-previous checkpoint if(fileExists(chkDataFileOld)){ if( ( std::remove( chkDataFileOld.c_str() ) + std::remove( chkParamFileOld.c_str() ) ) != 0 ) pcout << "WARNING: cannot remove old chekpoint file!" << std::endl; } // Rename previous checkpoint - if (std::rename(chkParamFile.c_str(), chkParamFileOld.c_str()) || std::rename(chkDataFile.c_str(), chkDataFileOld.c_str() )) + if (std::rename(chkParamFile.c_str(), chkParamFileOld.c_str()) || std::rename(chkDataFile.c_str(), chkDataFileOld.c_str() )) pcout << "WARNING: cannot rename old chekpoint file!" << std::endl; } } global::mpi().barrier(); - + plb_ofstream ofile(chkParamFile.c_str()); ofile << stat_cycle << endl; saveBinaryBlock(*lattice, chkDataFile); } diff --git a/io.cpp b/io.cpp index 5a7e873..a7e039b 100644 --- a/io.cpp +++ b/io.cpp @@ -18,7 +18,7 @@ int dirExists(const string& pathName) if( stat( pathName.c_str(), &info ) != 0 ) return -1; // Cannot acces path - else if( info.st_mode & S_IFDIR ) // S_ISDIR() doesn't exist on my windows + else if( info.st_mode & S_IFDIR ) // S_ISDIR() doesn't exist on my windows return 1; // Path exists else return 0; // Path does not exist @@ -94,7 +94,7 @@ void writeVTK(MultiBlockLattice3D& lattice, const SimPar &sim, pli vtkOut.writeData<6,float>(*computeShearStress(lattice), "sigma [1/m2s]", 1./(sim.C_l*sim.C_t*sim.C_t)); vtkOut.writeData(*computeSymmetricTensorNorm(*computeStrainRateFromStress(lattice)), "S_norm [1/s]", 1./sim.C_t ); // TODO - output viscosity? - + if (field1 != nullptr) vtkOut.writeData(*field1, "field1"); } @@ -114,7 +114,7 @@ void writeHDF5(MultiBlockLattice3D& lattice, const SimPar &sim, pl MultiScalarField3D DistributedDensity = *computeDensity(lattice); MultiTensorField3D DistributedShearStress = *computeShearStress(lattice); MultiScalarField3D DistributedS_Norm = *computeSymmetricTensorNorm(*computeStrainRateFromStress(lattice)); - // MultiScalarField3D DistributedField1 = *field1; // Used for additional fields, e.g. porosity, do any necessary calculations here. + // MultiScalarField3D DistributedField1 = *field1; // Used for additional fields, e.g. porosity, do any necessary calculations here. // Density/Velocity/... shared the same atomic block distribution! MultiBlockManagement3D VelocityBlockManagement = DistributedVelocity.getMultiBlockManagement(); @@ -145,7 +145,7 @@ void writeHDF5(MultiBlockLattice3D& lattice, const SimPar &sim, pl for(unsigned int i = LocalBulk.getBulk().x0; i <= LocalBulk.getBulk().x1; i++) for(unsigned int j = LocalBulk.getBulk().y0; j <= LocalBulk.getBulk().y1; j++) for(unsigned int k = LocalBulk.getBulk().z0; k <= LocalBulk.getBulk().z1; k++){ - + // Now we convert the global scale coordinates to block local coordinates unsigned int LocalX = LocalBulk.toLocalX(i); unsigned int LocalY = LocalBulk.toLocalY(j); @@ -173,11 +173,11 @@ void writeHDF5(MultiBlockLattice3D& lattice, const SimPar &sim, pl // Additional field - No unit conversion! if (field1 != nullptr) { double foundField1 = *field1->getComponent(blockId).get(LocalX, LocalY, LocalZ); - Field1.push_back(foundField1); + Field1.push_back(foundField1); } else { Field1.push_back(0.0); - } + } // Rank of current mpi thread this_rank.push_back(RankID); @@ -191,10 +191,10 @@ void writeHDF5(MultiBlockLattice3D& lattice, const SimPar &sim, pl ///////////////////////////// Saving HDF5 ///////////////////////////// - // Now save the partial local data to hdf5, if you dont understand, + // Now save the partial local data to hdf5, if you dont understand, // check (https://github.com/BlueBrain/HighFive/blob/master/src/examples/parallel_hdf5_collective_io.cpp) using namespace HighFive; - + FileAccessProps fapl; // Tell HDF5 to use MPI-IO fapl.add(MPIOFileAccess{MPI_COMM_WORLD, MPI_INFO_NULL}); @@ -208,8 +208,13 @@ void writeHDF5(MultiBlockLattice3D& lattice, const SimPar &sim, pl // For compression DataSetCreateProps props; - // Use chunking - props.add(Chunking(std::vector{100, 100, 100})); + // Use adaptive chunking based on dataset dimensions + // Chunk size should be <= dimension size + hsize_t chunk_x = std::min(Nx, 100); + hsize_t chunk_y = std::min(Ny, 100); + hsize_t chunk_z = std::min(Nz, 100); + // The order matters if chunking is not cubic + props.add(Chunking(std::vector{chunk_z, chunk_y, chunk_x})); // Enable shuffle props.add(Shuffle()); // Enable deflate @@ -233,7 +238,7 @@ void writeHDF5(MultiBlockLattice3D& lattice, const SimPar &sim, pl DataSet S_Norm = file.createDataSet("S_norm", DataSpace(Dims), props); // Field1 DataSet Field1_data = file.createDataSet("Field1", DataSpace(Dims), props); - + // MPI rank DataSet Rank = file.createDataSet("MPI_rank", DataSpace(Dims), props); @@ -416,7 +421,7 @@ void writeNPZ(MultiBlockLattice3D& lattice, plint iter) if(global::mpi().isMainProcessor()) { double *data = new double[3*nx*ny*nz]; - for(unsigned int i = 0; i < nx; i++) + for(unsigned int i = 0; i < nx; i++) for(unsigned int j = 0; j < ny; j++) for(unsigned int k = 0; k < nz; k++) { int idx = (i*nx*nz+j*nz+k)*3; @@ -426,7 +431,7 @@ void writeNPZ(MultiBlockLattice3D& lattice, plint iter) data[idx+2] = localVelocity.get(i, j, k)[2]; } - cnpy::npz_save(createFileName("output_", iter, 6) + ".npz", "velocity",&data[0],{3,nz,ny,nx},"w"); + cnpy::npz_save(createFileName("output_", iter, 6) + ".npz", "velocity",&data[0],{3,nz,ny,nx},"w"); } global::mpi().barrier(); diff --git a/testcase/input/simple_flowrate.txt b/testcase/input/simple_flowrate.txt new file mode 100644 index 0000000..23ccf3d --- /dev/null +++ b/testcase/input/simple_flowrate.txt @@ -0,0 +1,286 @@ +285 +0.000000 0.300000 +0.020000 0.377778 +0.040000 0.455556 +0.060000 0.533333 +0.080000 0.611111 +0.100000 0.688889 +0.120000 0.766667 +0.140000 0.844444 +0.160000 0.922222 +0.180000 1.000000 +0.200000 1.000000 +0.220000 1.000000 +0.240000 1.000000 +0.260000 1.000000 +0.280000 1.000000 +0.300000 1.000000 +0.320000 0.979412 +0.340000 0.958824 +0.360000 0.938235 +0.380000 0.917647 +0.400000 0.897059 +0.420000 0.876471 +0.440000 0.855882 +0.460000 0.835294 +0.480000 0.814706 +0.500000 0.794118 +0.520000 0.773529 +0.540000 0.752941 +0.560000 0.732353 +0.580000 0.711765 +0.600000 0.691176 +0.620000 0.670588 +0.640000 0.650000 +0.660000 0.629412 +0.680000 0.608824 +0.700000 0.588235 +0.720000 0.567647 +0.740000 0.547059 +0.760000 0.526471 +0.780000 0.505882 +0.800000 0.485294 +0.820000 0.464706 +0.840000 0.444118 +0.860000 0.423529 +0.880000 0.402941 +0.900000 0.382353 +0.920000 0.361765 +0.940000 0.341176 +0.960000 0.320588 +0.980000 0.300000 +1.000000 0.300000 +1.020000 0.500000 +1.040000 0.700000 +1.060000 0.700000 +1.080000 0.700000 +1.100000 0.660000 +1.120000 0.620000 +1.140000 0.580000 +1.160000 0.540000 +1.180000 0.500000 +1.200000 0.460000 +1.220000 0.420000 +1.240000 0.380000 +1.260000 0.340000 +1.280000 0.300000 +1.300000 0.300000 +1.320000 0.500000 +1.340000 0.700000 +1.360000 0.700000 +1.380000 0.700000 +1.400000 0.660000 +1.420000 0.620000 +1.440000 0.580000 +1.460000 0.540000 +1.480000 0.500000 +1.500000 0.460000 +1.520000 0.420000 +1.540000 0.380000 +1.560000 0.340000 +1.580000 0.300000 +1.600000 0.300000 +1.620000 0.500000 +1.640000 0.700000 +1.660000 0.700000 +1.680000 0.700000 +1.700000 0.660000 +1.720000 0.620000 +1.740000 0.580000 +1.760000 0.540000 +1.780000 0.500000 +1.800000 0.460000 +1.820000 0.420000 +1.840000 0.380000 +1.860000 0.340000 +1.880000 0.300000 +1.900000 0.300000 +1.920000 0.377778 +1.940000 0.455556 +1.960000 0.533333 +1.980000 0.611111 +2.000000 0.688889 +2.020000 0.766667 +2.040000 0.844444 +2.060000 0.922222 +2.080000 1.000000 +2.100000 1.000000 +2.120000 1.000000 +2.140000 1.000000 +2.160000 1.000000 +2.180000 1.000000 +2.200000 1.000000 +2.220000 0.979412 +2.240000 0.958824 +2.260000 0.938235 +2.280000 0.917647 +2.300000 0.897059 +2.320000 0.876471 +2.340000 0.855882 +2.360000 0.835294 +2.380000 0.814706 +2.400000 0.794118 +2.420000 0.773529 +2.440000 0.752941 +2.460000 0.732353 +2.480000 0.711765 +2.500000 0.691176 +2.520000 0.670588 +2.540000 0.650000 +2.560000 0.629412 +2.580000 0.608824 +2.600000 0.588235 +2.620000 0.567647 +2.640000 0.547059 +2.660000 0.526471 +2.680000 0.505882 +2.700000 0.485294 +2.720000 0.464706 +2.740000 0.444118 +2.760000 0.423529 +2.780000 0.402941 +2.800000 0.382353 +2.820000 0.361765 +2.840000 0.341176 +2.860000 0.320588 +2.880000 0.300000 +2.900000 0.300000 +2.920000 0.500000 +2.940000 0.700000 +2.960000 0.700000 +2.980000 0.700000 +3.000000 0.660000 +3.020000 0.620000 +3.040000 0.580000 +3.060000 0.540000 +3.080000 0.500000 +3.100000 0.460000 +3.120000 0.420000 +3.140000 0.380000 +3.160000 0.340000 +3.180000 0.300000 +3.200000 0.300000 +3.220000 0.500000 +3.240000 0.700000 +3.260000 0.700000 +3.280000 0.700000 +3.300000 0.660000 +3.320000 0.620000 +3.340000 0.580000 +3.360000 0.540000 +3.380000 0.500000 +3.400000 0.460000 +3.420000 0.420000 +3.440000 0.380000 +3.460000 0.340000 +3.480000 0.300000 +3.500000 0.300000 +3.520000 0.500000 +3.540000 0.700000 +3.560000 0.700000 +3.580000 0.700000 +3.600000 0.660000 +3.620000 0.620000 +3.640000 0.580000 +3.660000 0.540000 +3.680000 0.500000 +3.700000 0.460000 +3.720000 0.420000 +3.740000 0.380000 +3.760000 0.340000 +3.780000 0.300000 +3.800000 0.300000 +3.820000 0.377778 +3.840000 0.455556 +3.860000 0.533333 +3.880000 0.611111 +3.900000 0.688889 +3.920000 0.766667 +3.940000 0.844444 +3.960000 0.922222 +3.980000 1.000000 +4.000000 1.000000 +4.020000 1.000000 +4.040000 1.000000 +4.060000 1.000000 +4.080000 1.000000 +4.100000 1.000000 +4.120000 0.979412 +4.140000 0.958824 +4.160000 0.938235 +4.180000 0.917647 +4.200000 0.897059 +4.220000 0.876471 +4.240000 0.855882 +4.260000 0.835294 +4.280000 0.814706 +4.300000 0.794118 +4.320000 0.773529 +4.340000 0.752941 +4.360000 0.732353 +4.380000 0.711765 +4.400000 0.691176 +4.420000 0.670588 +4.440000 0.650000 +4.460000 0.629412 +4.480000 0.608824 +4.500000 0.588235 +4.520000 0.567647 +4.540000 0.547059 +4.560000 0.526471 +4.580000 0.505882 +4.600000 0.485294 +4.620000 0.464706 +4.640000 0.444118 +4.660000 0.423529 +4.680000 0.402941 +4.700000 0.382353 +4.720000 0.361765 +4.740000 0.341176 +4.760000 0.320588 +4.780000 0.300000 +4.800000 0.300000 +4.820000 0.500000 +4.840000 0.700000 +4.860000 0.700000 +4.880000 0.700000 +4.900000 0.660000 +4.920000 0.620000 +4.940000 0.580000 +4.960000 0.540000 +4.980000 0.500000 +5.000000 0.460000 +5.020000 0.420000 +5.040000 0.380000 +5.060000 0.340000 +5.080000 0.300000 +5.100000 0.300000 +5.120000 0.500000 +5.140000 0.700000 +5.160000 0.700000 +5.180000 0.700000 +5.200000 0.660000 +5.220000 0.620000 +5.240000 0.580000 +5.260000 0.540000 +5.280000 0.500000 +5.300000 0.460000 +5.320000 0.420000 +5.340000 0.380000 +5.360000 0.340000 +5.380000 0.300000 +5.400000 0.300000 +5.420000 0.500000 +5.440000 0.700000 +5.460000 0.700000 +5.480000 0.700000 +5.500000 0.660000 +5.520000 0.620000 +5.540000 0.580000 +5.560000 0.540000 +5.580000 0.500000 +5.600000 0.460000 +5.620000 0.420000 +5.640000 0.380000 +5.660000 0.340000 +5.680000 0.300000 diff --git a/testcase/input/simple_tube_small.npz b/testcase/input/simple_tube_small.npz new file mode 100644 index 0000000..95db4e7 Binary files /dev/null and b/testcase/input/simple_tube_small.npz differ diff --git a/testcase/test_SimpleTube_Freeflow.xml b/testcase/test_SimpleTube_Freeflow.xml new file mode 100644 index 0000000..c1ce0cd --- /dev/null +++ b/testcase/test_SimpleTube_Freeflow.xml @@ -0,0 +1,48 @@ + +
1e-4
+ 30 + output_simple_tube + + 3.2 + 0.016 + + + + + + ./input/simple_tube_small.npz + + + + + Inlet + 1 + + ./input/simple_flowrate.txt + + 0.1 + + + + Outlet 0 Freeflow + 4 + + + + + + + + Outlet 2 Pressure + 3 + + + 0 + + + + + + 100 + 200 +