Skip to content

Commit

Permalink
Merge pull request #145 from AMReX-FHD/disjoining
Browse files Browse the repository at this point in the history
Add disjoining pressure
  • Loading branch information
ajnonaka authored Sep 12, 2024
2 parents 19dfab0 + f1d65bf commit 5a2ef5c
Show file tree
Hide file tree
Showing 3 changed files with 47 additions and 4 deletions.
47 changes: 43 additions & 4 deletions exec/thinFilm/main_driver.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -137,7 +137,9 @@ void main_driver(const char* argv)

MultiFab height(ba, dmap, 1, 1);
MultiFab Laph (ba, dmap, 1, 1);
MultiFab disjoining(ba, dmap, 1, 1);
Laph.setVal(0.); // prevent intermediate NaN calculations behind physical boundaries
disjoining.setVal(0.);

// for statsitics

Expand Down Expand Up @@ -169,6 +171,11 @@ void main_driver(const char* argv)
gradLaph[1].define(convert(ba,nodal_flag_y), dmap, 1, 0);,
gradLaph[2].define(convert(ba,nodal_flag_z), dmap, 1, 0););

std::array< MultiFab, AMREX_SPACEDIM > gradDisjoining;
AMREX_D_TERM(gradDisjoining[0].define(convert(ba,nodal_flag_x), dmap, 1, 0);,
gradDisjoining[1].define(convert(ba,nodal_flag_y), dmap, 1, 0);,
gradDisjoining[2].define(convert(ba,nodal_flag_z), dmap, 1, 0););

std::array< MultiFab, AMREX_SPACEDIM > flux;
AMREX_D_TERM(flux[0] .define(convert(ba,nodal_flag_x), dmap, 1, 0);,
flux[1] .define(convert(ba,nodal_flag_y), dmap, 1, 0);,
Expand Down Expand Up @@ -205,6 +212,7 @@ void main_driver(const char* argv)
// constant factor in noise term
Real ConstNoise = 2.*k_B*T_init[0] / (3.*visc_coef);
Real Const3dx = thinfilm_gamma / (3.*visc_coef);
Real Const3dx_nogamma = 1. / (3.*visc_coef);

Real time = 0.;

Expand Down Expand Up @@ -365,12 +373,16 @@ void main_driver(const char* argv)

}

// compute Laph
// compute Laph and disjoining
for ( MFIter mfi(Laph,TilingIfNotGPU()); mfi.isValid(); ++mfi ) {

const Box& bx = mfi.tilebox();

const Array4<Real> & L = Laph.array(mfi);

const Array4<Real> & h = height.array(mfi);

const Array4<Real> & Disjoining = disjoining.array(mfi);

AMREX_D_TERM(const Array4<Real> & gradhx = gradh[0].array(mfi);,
const Array4<Real> & gradhy = gradh[1].array(mfi);,
Expand All @@ -380,11 +392,14 @@ void main_driver(const char* argv)
{
L(i,j,k) = x_flux_fac * (gradhx(i+1,j,k) - gradhx(i,j,k)) / dx[0]
+ y_flux_fac * (gradhy(i,j+1,k) - gradhy(i,j,k)) / dx[1];

Disjoining(i,j,k) = thinfilm_hamaker / (6.*M_PI*std::pow(h(i,j,k),3.));
});
}
Laph.FillBoundary(geom.periodicity());
disjoining.FillBoundary(geom.periodicity());

// compute gradLaph
// compute gradLaph and gradDisjoining
for ( MFIter mfi(height,TilingIfNotGPU()); mfi.isValid(); ++mfi ) {

AMREX_D_TERM(const Box & bx_x = mfi.nodaltilebox(0);,
Expand All @@ -394,17 +409,25 @@ void main_driver(const char* argv)
AMREX_D_TERM(const Array4<Real> & gradLaphx = gradLaph[0].array(mfi);,
const Array4<Real> & gradLaphy = gradLaph[1].array(mfi);,
const Array4<Real> & gradLaphz = gradLaph[2].array(mfi););

AMREX_D_TERM(const Array4<Real> & gradDisjoiningx = gradDisjoining[0].array(mfi);,
const Array4<Real> & gradDisjoiningy = gradDisjoining[1].array(mfi);,
const Array4<Real> & gradDisjoiningz = gradDisjoining[2].array(mfi););

const Array4<Real> & L = Laph.array(mfi);

const Array4<Real> & Disjoining = disjoining.array(mfi);

amrex::ParallelFor(bx_x, bx_y,
[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
gradLaphx(i,j,k) = ( L(i,j,k) - L(i-1,j,k) ) / dx[0];
gradDisjoiningx(i,j,k) = ( Disjoining(i,j,k) -Disjoining(i-1,j,k) ) / dx[0];
},
[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
gradLaphy(i,j,k) = ( L(i,j,k) - L(i,j-1,k) ) / dx[1];
gradDisjoiningy(i,j,k) = ( Disjoining(i,j,k) -Disjoining(i,j-1,k) ) / dx[1];
});
}

Expand All @@ -431,18 +454,24 @@ void main_driver(const char* argv)
const Array4<Real> & randfacey = randface[1].array(mfi);,
const Array4<Real> & randfacez = randface[2].array(mfi););

AMREX_D_TERM(const Array4<Real> & gradDisjoiningx = gradDisjoining[0].array(mfi);,
const Array4<Real> & gradDisjoiningy = gradDisjoining[1].array(mfi);,
const Array4<Real> & gradDisjoiningz = gradDisjoining[2].array(mfi););

amrex::ParallelFor(bx_x, bx_y,
[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
fluxx(i,j,k) = x_flux_fac * (
std::sqrt(ConstNoise*std::pow(hfacex(i,j,k),3.) / (dt*dVol)) * randfacex(i,j,k)
+ Const3dx * std::pow(hfacex(i,j,k),3.)*gradLaphx(i,j,k) );
+ Const3dx * std::pow(hfacex(i,j,k),3.)*gradLaphx(i,j,k)
+ Const3dx_nogamma * std::pow(hfacex(i,j,k),3.)*gradDisjoiningx(i,j,k));
},
[=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{
fluxy(i,j,k) = y_flux_fac * (
std::sqrt(ConstNoise*std::pow(hfacey(i,j,k),3.) / (dt*dVol)) * randfacey(i,j,k)
+ Const3dx * std::pow(hfacey(i,j,k),3.)*gradLaphy(i,j,k) );
+ Const3dx * std::pow(hfacey(i,j,k),3.)*gradLaphy(i,j,k)
+ Const3dx_nogamma * std::pow(hfacex(i,j,k),3.)*gradDisjoiningy(i,j,k) );
});

// lo x-faces
Expand Down Expand Up @@ -489,8 +518,11 @@ void main_driver(const char* argv)

const Array4<Real> & h = height.array(mfi);

// amrex::Print{} << "HEIGHT " << time << " " << h(0,0,0) << " " << h(31,0,0) << std::endl;

amrex::ParallelFor(bx, [=] AMREX_GPU_DEVICE (int i, int j, int k) noexcept
{

h(i,j,k) -= dt * ( (fluxx(i+1,j,k) - fluxx(i,j,k)) / dx[0]
+(fluxy(i,j+1,k) - fluxy(i,j,k)) / dx[1]);
});
Expand All @@ -511,6 +543,8 @@ void main_driver(const char* argv)
// copy distributed data into 1D data
height_onegrid.ParallelCopy(height, 0, 0, 1);

amrex::Real sumh = 0.;

for ( MFIter mfi(height_onegrid,false); mfi.isValid(); ++mfi ) {

std::ofstream hstream;
Expand All @@ -530,9 +564,14 @@ void main_driver(const char* argv)
for (auto j = lo.y; j <= hi.y; ++j) {
for (auto i = lo.x; i <= hi.x; ++i) {
hstream << std::setprecision(15) << mfdata(i,j,0) << " ";
sumh += mfdata(i,j,0);
if(j==0 && i==0){
amrex::Print{} << "HEIGHT " << time << " " << mfdata(0,0,0) << " " << mfdata(31,0,0) << std::endl;
}
}
hstream << "\n";
}
amrex::Print{} << "SUM " << sumh << std::endl;

} // end MFIter

Expand Down
3 changes: 3 additions & 0 deletions exec/thinFilm/thinfilm_functions.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,7 @@
AMREX_GPU_MANAGED amrex::Real thinfilm::thinfilm_h0;
AMREX_GPU_MANAGED amrex::Real thinfilm::thinfilm_gamma;
AMREX_GPU_MANAGED amrex::Real thinfilm::thinfilm_pertamp;
AMREX_GPU_MANAGED amrex::Real thinfilm::thinfilm_hamaker;

AMREX_GPU_MANAGED int thinfilm::thinfilm_icorr;
AMREX_GPU_MANAGED int thinfilm::thinfilm_jcorr;
Expand All @@ -19,6 +20,7 @@ void InitializeThinfilmNamespace() {
thinfilm_icorr = 0;
thinfilm_jcorr = 0;
thinfilm_pertamp = 0.;
thinfilm_hamaker = 0.;

do_fft_diag = 1;

Expand All @@ -28,6 +30,7 @@ void InitializeThinfilmNamespace() {
pp.get("thinfilm_gamma",thinfilm_gamma);

pp.query("thinfilm_pertamp",thinfilm_pertamp);
pp.query("thinfilm_hamaker",thinfilm_hamaker);

pp.query("thinfilm_icorr",thinfilm_icorr);
pp.query("thinfilm_jcorr",thinfilm_jcorr);
Expand Down
1 change: 1 addition & 0 deletions exec/thinFilm/thinfilm_namespace.H
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@ namespace thinfilm {
extern AMREX_GPU_MANAGED amrex::Real thinfilm_h0;
extern AMREX_GPU_MANAGED amrex::Real thinfilm_gamma;
extern AMREX_GPU_MANAGED amrex::Real thinfilm_pertamp;
extern AMREX_GPU_MANAGED amrex::Real thinfilm_hamaker;

extern AMREX_GPU_MANAGED int thinfilm_icorr;
extern AMREX_GPU_MANAGED int thinfilm_jcorr;
Expand Down

0 comments on commit 5a2ef5c

Please sign in to comment.