Thanks for your help in advance,

I am trying to implement the Cahn-Hilliard equation in PALABOS to simulate

multiphase flow with high density ratios (see Banari et. al implementation for

reference “An improved two-phase Lattice Boltzmann model for high density

ratios: application to wave breaking”). The equilibrium function is expressed as

(equation 21):

The first and final terms are not problematic, but the second term has a second

order spatial derivative with respect to $\phi$ (equation 2 in Banari et al.):

From my understanding all **dynamics** objects have a *computeEquilibrium*

method that operates cell-wise and always takes these arguments (from

src/core/dynamics.h):

```
virtual void computeEquilibrium(plint iPop, T rhoBar,
Array<T,Descriptor<T>::d> const& j, T jSqr, T thetaBar=T())
```

My approach was to have an additional term in this expression (e.g. *T lap_rho*)

that passes the laplacian as an argument to *computeEquilibrium*, but I am

unable to find where the *computeEquilibrium* method is actually called.

Put another way, I am trying to modify the call to *computeEquilibirium* so that

it can have the spatial derivatives as another argument.

I assume this call is integrated into the *collideAndStream* method of the

MultiBlockLattice, but extending the MultiBlockLattice seems to be going against

the recommendations of the User Guide.

Thanks again for help.

Cyrus