From f02cf7abc491c7a4bb6faecb33d4d58b5dbd1d78 Mon Sep 17 00:00:00 2001 From: Rex Zhe Li Date: Sun, 28 Mar 2021 21:26:01 -0400 Subject: [PATCH 1/2] remove the phase field denoise since it ruins mass conservation --- cpu/FreeLee.cpp | 6 ------ 1 file changed, 6 deletions(-) diff --git a/cpu/FreeLee.cpp b/cpu/FreeLee.cpp index 1e2fc04d..5daac030 100644 --- a/cpu/FreeLee.cpp +++ b/cpu/FreeLee.cpp @@ -161,9 +161,6 @@ extern "C" void ScaLBL_D3Q7_AAodd_FreeLeeModel_PhaseField(int *neighborList, int fq = hq[nread]; phi += fq; - if (phi > 1.f) phi = 1.0; - if (phi < -1.f) phi = -1.0; - // save the number densities Den[n] = rhoA + 0.5*(1.0-phi)*(rhoB-rhoA); @@ -210,9 +207,6 @@ extern "C" void ScaLBL_D3Q7_AAeven_FreeLeeModel_PhaseField(int *Map, double *hq, fq = hq[5*Np+n]; phi += fq; - if (phi > 1.f) phi = 1.0; - if (phi < -1.f) phi = -1.0; - // save the number densities Den[n] = rhoA + 0.5*(1.0-phi)*(rhoB-rhoA); From e87db53a019bb93566020e035950e0fbb8e2fb56 Mon Sep 17 00:00:00 2001 From: Rex Zhe Li Date: Mon, 29 Mar 2021 23:28:26 -0400 Subject: [PATCH 2/2] apply the phase denoise without affecting the mass conservation --- cpu/FreeLee.cpp | 30 ++++++++++++++++++++++-------- 1 file changed, 22 insertions(+), 8 deletions(-) diff --git a/cpu/FreeLee.cpp b/cpu/FreeLee.cpp index 5daac030..98613f5c 100644 --- a/cpu/FreeLee.cpp +++ b/cpu/FreeLee.cpp @@ -1523,6 +1523,7 @@ extern "C" void ScaLBL_D3Q19_AAodd_FreeLeeModel_Combined(int *neighborList, int double tau;//position dependent LB relaxation time for fluid double C,theta; double M = 2.0/9.0*(tauM-0.5);//diffusivity (or mobility) for the phase field D3Q7 + double phi_temp; for (int n=start; n1.f) phi_temp=1.0; + if (phi<-1.f) phi_temp=-1.0; // local relaxation time tau=tauA + 0.5*(1.0-phi)*(tauB-tauA); @@ -1673,8 +1677,10 @@ extern "C" void ScaLBL_D3Q19_AAodd_FreeLeeModel_Combined(int *neighborList, int ny = -3.0*1.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18)); nz = -3.0*1.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18)); //............Compute the Chemical Potential............................... - chem = 2.0*3.0/18.0*(m1+m2+m3+m4+m5+m6-6*phi+0.5*(m7+m8+m9+m10+m11+m12+m13+m14+m15+m16+m17+m18-12*phi));//intermediate var, i.e. the laplacian - chem = 4.0*beta*phi*(phi+1.0)*(phi-1.0)-kappa*chem; + //chem = 2.0*3.0/18.0*(m1+m2+m3+m4+m5+m6-6*phi+0.5*(m7+m8+m9+m10+m11+m12+m13+m14+m15+m16+m17+m18-12*phi));//intermediate var, i.e. the laplacian + //chem = 4.0*beta*phi*(phi+1.0)*(phi-1.0)-kappa*chem; + chem = 2.0*3.0/18.0*(m1+m2+m3+m4+m5+m6-6*phi_temp+0.5*(m7+m8+m9+m10+m11+m12+m13+m14+m15+m16+m17+m18-12*phi_temp));//intermediate var, i.e. the laplacian + chem = 4.0*beta*phi_temp*(phi_temp+1.0)*(phi_temp-1.0)-kappa*chem; //............Compute the Mixed Gradient................................... mgx = -3.0*1.0/18.0*(mm1-mm2+0.5*(mm7-mm8+mm9-mm10+mm11-mm12+mm13-mm14)); mgy = -3.0*1.0/18.0*(mm3-mm4+0.5*(mm7-mm8-mm9+mm10+mm15-mm16+mm17-mm18)); @@ -1685,7 +1691,7 @@ extern "C" void ScaLBL_D3Q19_AAodd_FreeLeeModel_Combined(int *neighborList, int if (C<1.0e-12) nx=ny=nz=0.0; double mg_mag = sqrt(mgx*mgx+mgy*mgy+mgz*mgz); if (mg_mag<1.0e-12) mgx=mgy=mgz=0.0; - //TODO - maybe you can also de-noise chemical potential ? within the bulk phase chem should be ZERO + //maybe you can also de-noise chemical potential ? within the bulk phase chem should be ZERO if (fabs(chem)<1.0e-12) chem=0.0; // q=0 @@ -2052,7 +2058,8 @@ extern "C" void ScaLBL_D3Q19_AAodd_FreeLeeModel_Combined(int *neighborList, int ny = ny/ColorMag; nz = nz/ColorMag; //compute surface tension-related parameter - theta = 4.5*M*2.0*(1-phi*phi)/W; + //theta = 4.5*M*2.0*(1-phi*phi)/W; + theta = 4.5*M*2.0*(1-phi_temp*phi_temp)/W; //load distributions of phase field //q=0 @@ -2139,6 +2146,7 @@ extern "C" void ScaLBL_D3Q19_AAeven_FreeLeeModel_Combined(int *Map, double *dist double tau;//position dependent LB relaxation time for fluid double C,theta; double M = 2.0/9.0*(tauM-0.5);//diffusivity (or mobility) for the phase field D3Q7 + double phi_temp; for (int n=start; n1.f) phi_temp=1.0; + if (phi<-1.f) phi_temp=-1.0; // local relaxation time tau=tauA + 0.5*(1.0-phi)*(tauB-tauA); @@ -2289,8 +2300,10 @@ extern "C" void ScaLBL_D3Q19_AAeven_FreeLeeModel_Combined(int *Map, double *dist ny = -3.0*1.0/18.0*(m3-m4+0.5*(m7-m8-m9+m10+m15-m16+m17-m18)); nz = -3.0*1.0/18.0*(m5-m6+0.5*(m11-m12-m13+m14+m15-m16-m17+m18)); //............Compute the Chemical Potential............................... - chem = 2.0*3.0/18.0*(m1+m2+m3+m4+m5+m6-6*phi+0.5*(m7+m8+m9+m10+m11+m12+m13+m14+m15+m16+m17+m18-12*phi));//intermediate var, i.e. the laplacian - chem = 4.0*beta*phi*(phi+1.0)*(phi-1.0)-kappa*chem; + //chem = 2.0*3.0/18.0*(m1+m2+m3+m4+m5+m6-6*phi+0.5*(m7+m8+m9+m10+m11+m12+m13+m14+m15+m16+m17+m18-12*phi));//intermediate var, i.e. the laplacian + //chem = 4.0*beta*phi*(phi+1.0)*(phi-1.0)-kappa*chem; + chem = 2.0*3.0/18.0*(m1+m2+m3+m4+m5+m6-6*phi_temp+0.5*(m7+m8+m9+m10+m11+m12+m13+m14+m15+m16+m17+m18-12*phi_temp));//intermediate var, i.e. the laplacian + chem = 4.0*beta*phi_temp*(phi_temp+1.0)*(phi_temp-1.0)-kappa*chem; //............Compute the Mixed Gradient................................... mgx = -3.0*1.0/18.0*(mm1-mm2+0.5*(mm7-mm8+mm9-mm10+mm11-mm12+mm13-mm14)); mgy = -3.0*1.0/18.0*(mm3-mm4+0.5*(mm7-mm8-mm9+mm10+mm15-mm16+mm17-mm18)); @@ -2301,7 +2314,7 @@ extern "C" void ScaLBL_D3Q19_AAeven_FreeLeeModel_Combined(int *Map, double *dist if (C<1.0e-12) nx=ny=nz=0.0; double mg_mag = sqrt(mgx*mgx+mgy*mgy+mgz*mgz); if (mg_mag<1.0e-12) mgx=mgy=mgz=0.0; - //TODO - maybe you can also de-noise chemical potential ? within the bulk phase chem should be ZERO + //maybe you can also de-noise chemical potential ? within the bulk phase chem should be ZERO if (fabs(chem)<1.0e-12) chem=0.0; // q=0 @@ -2651,7 +2664,8 @@ extern "C" void ScaLBL_D3Q19_AAeven_FreeLeeModel_Combined(int *Map, double *dist ny = ny/ColorMag; nz = nz/ColorMag; //compute surface tension-related parameter - theta = 4.5*M*2.0*(1-phi*phi)/W; + //theta = 4.5*M*2.0*(1-phi*phi)/W; + theta = 4.5*M*2.0*(1-phi_temp*phi_temp)/W; //load distributions of phase field //q=0