OilfieldToolsNet/LBPM
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

merge master with crusher

Browse Source
This commit is contained in:
James E McClure 2022-03-27 16:02:56 -04:00
commit 6911e97c29
7 changed files with 30 additions and 12 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
analysis/Minkowski.cpp
View File

@ -67,6 +67,7 @@ Minkowski::~Minkowski() {
void Minkowski::ComputeScalar(const DoubleArray &Field, const double isovalue) {
PROFILE_START("ComputeScalar");
Xi = Ji = Ai = 0.0;
DCEL object;
int e1, e2, e3;
@ -160,6 +161,7 @@ void Minkowski::MeasureObject() {
* 1 - labels the rest of the
*/
//DoubleArray smooth_distance(Nx,Ny,Nz);
for (int k = 0; k < Nz; k++) {
for (int j = 0; j < Ny; j++) {
for (int i = 0; i < Nx; i++) {
@ -168,6 +170,8 @@ void Minkowski::MeasureObject() {
}
}
CalcDist(distance, id, *Dm);
Dm->CommunicateMeshHalo(distance);
//Mean3D(distance,smooth_distance);
//Eikonal(distance, id, *Dm, 20, {true, true, true});
ComputeScalar(distance, 0.0);
@ -179,7 +183,7 @@ void Minkowski::MeasureObject(double factor, const DoubleArray &Phi) {
*
* THIS ALGORITHM ASSUMES THAT id() is populated with phase id to distinguish objects
* 0 - labels the object
* 1 - labels the rest of the
* 1 - labels the rest
*/
for (int k = 0; k < Nz; k++) {
for (int j = 0; j < Ny; j++) {

8
analysis/SubPhase.cpp
View File

@ -595,7 +595,7 @@ void SubPhase::Full() {
for (j = 0; j < Ny; j++) {
for (i = 0; i < Nx; i++) {
n = k * Nx * Ny + j * Nx + i;
if (!(Dm->id[n] > 0)) {
if (SDs(n) <= 0.0) {
// Solid phase
morph_n->id(i, j, k) = 1;
@ -642,7 +642,7 @@ void SubPhase::Full() {
for (j = 0; j < Ny; j++) {
for (i = 0; i < Nx; i++) {
n = k * Nx * Ny + j * Nx + i;
if (!(Dm->id[n] > 0)) {
if (SDs(n) <= 0.0) {
// Solid phase
morph_w->id(i, j, k) = 1;
@ -688,7 +688,7 @@ void SubPhase::Full() {
for (j = 0; j < Ny; j++) {
for (i = 0; i < Nx; i++) {
n = k * Nx * Ny + j * Nx + i;
if (!(Dm->id[n] > 0)) {
if (SDs(n) <= 0.0) {
// Solid phase
morph_i->id(i, j, k) = 1;
} else if (DelPhi(n) > 1e-4) {
@ -731,7 +731,7 @@ void SubPhase::Full() {
for (i = imin; i < Nx - 1; i++) {
n = k * Nx * Ny + j * Nx + i;
// Compute volume averages
if (Dm->id[n] > 0) {
if (SDs(n) > 0.0) {
// compute density
double nA = Rho_n(n);
double nB = Rho_w(n);

4
models/IonModel.cpp
View File

@ -170,7 +170,7 @@ void ScaLBL_IonModel::ReadParams(string filename, vector<int> &num_iter) {
BoundaryConditionSolid = ion_db->getScalar<int>("BC_Solid");
}
// Read boundary condition for ion transport
// BC = 0: normal periodic BC
// BC = 0: zero-flux bounce-back BC
// BC = 1: fixed ion concentration; unit=[mol/m^3]
// BC = 2: fixed ion flux (inward flux); unit=[mol/m^2/sec]
BoundaryConditionInlet.push_back(0);
@ -428,7 +428,7 @@ void ScaLBL_IonModel::ReadParams(string filename) {
BoundaryConditionSolid = ion_db->getScalar<int>("BC_Solid");
}
// Read boundary condition for ion transport
// BC = 0: normal periodic BC
// BC = 0: zero-flux bounce-back BC
// BC = 1: fixed ion concentration; unit=[mol/m^3]
// BC = 2: fixed ion flux (inward flux); unit=[mol/m^2/sec]
BoundaryConditionInlet.push_back(0);

3
models/MultiPhysController.cpp
View File

@ -4,7 +4,7 @@ ScaLBL_Multiphys_Controller::ScaLBL_Multiphys_Controller(
int RANK, int NP, const Utilities::MPI &COMM)
: rank(RANK), nprocs(NP), Restart(0), timestepMax(0), num_iter_Stokes(0),
num_iter_Ion(0), analysis_interval(0), visualization_interval(0),
tolerance(0), time_conv_max(0), comm(COMM) {}
tolerance(0), time_conv_max(0), time_conv_MainLoop(0), comm(COMM) {}
ScaLBL_Multiphys_Controller::~ScaLBL_Multiphys_Controller() {}
void ScaLBL_Multiphys_Controller::ReadParams(string filename) {
@ -22,6 +22,7 @@ void ScaLBL_Multiphys_Controller::ReadParams(string filename) {
visualization_interval = 10000;
tolerance = 1.0e-6;
time_conv_max = 0.0;
time_conv_MainLoop = 0.0;
// load input parameters
if (study_db->keyExists("timestepMax")) {

1
models/MultiPhysController.h
View File

@ -40,6 +40,7 @@ public:
int visualization_interval;
double tolerance;
double time_conv_max;
double time_conv_MainLoop;
//double SchmidtNum;//Schmidt number = kinematic_viscosity/mass_diffusivity
int rank, nprocs;

1
tests/TestTorus.cpp
View File

@ -168,6 +168,7 @@ int main(int argc, char **argv)
// Averages->Reduce();
Object->MeasureObject();
//Object->MeasureConnectedPathway();
double Vi = Object->V();
double Ai = Object->A();
double Hi = Object->H();

19
tests/lbpm_electrokinetic_SingleFluid_simulator.cpp
View File

@ -79,17 +79,29 @@ int main(int argc, char **argv)
IonModel.timestepMax = Study.getIonNumIter_PNP_coupling(StokesModel.time_conv,IonModel.time_conv);
IonModel.Initialize();
// Get maximal time converting factor based on Sotkes and Ion solvers
Study.getTimeConvMax_PNP_coupling(StokesModel.time_conv,IonModel.time_conv);
//Study.getTimeConvMax_PNP_coupling(StokesModel.time_conv,IonModel.time_conv);
// Get time conversion factor for the main iteration loop in electrokinetic single fluid simulator
Study.time_conv_MainLoop = StokesModel.timestepMax*StokesModel.time_conv;
// Initialize LB-Poisson model
PoissonSolver.ReadParams(filename);
PoissonSolver.SetDomain();
PoissonSolver.ReadInput();
PoissonSolver.Create();
PoissonSolver.Initialize(Study.time_conv_max);
PoissonSolver.Initialize(Study.time_conv_MainLoop);
if (rank == 0){
printf("********************************************************\n");
printf("Key Summary of LBPM electrokinetic single-fluid solver \n");
printf(" 1. Max LB Timestep: %i [lt]\n", Study.timestepMax);
printf(" 2. Time conversion factor per LB Timestep: %.6g [sec/lt]\n",Study.time_conv_MainLoop);
printf(" 3. Max Physical Time: %.6g [sec]\n",Study.timestepMax*Study.time_conv_MainLoop);
printf("********************************************************\n");
}
int timestep=0;
while (timestep < Study.timestepMax){
@ -98,7 +110,6 @@ int main(int argc, char **argv)
StokesModel.Run_Lite(IonModel.ChargeDensity, PoissonSolver.ElectricField);// Solve the N-S equations to get velocity
IonModel.Run(StokesModel.Velocity,PoissonSolver.ElectricField); //solve for ion transport and electric potential
timestep++;//AA operations
if (timestep%Study.analysis_interval==0){
Analysis.Basic(IonModel,PoissonSolver,StokesModel,timestep);
@ -116,7 +127,7 @@ int main(int argc, char **argv)
if (rank==0) printf("Save simulation raw data at maximum timestep\n");
Analysis.WriteVis(IonModel,PoissonSolver,StokesModel,Study.db,timestep);
if (rank==0) printf("Maximum timestep is reached and the simulation is completed\n");
if (rank==0) printf("Maximum LB timestep = %i is reached and the simulation is completed\n",Study.timestepMax);
if (rank==0) printf("*************************************************************\n");
PROFILE_STOP("Main");