[amwitt]

Hi all,

I am looking to add a variable diffusion coefficient to an advection diffusion equation I added in interFoam to model mass transport in turbulent free surface flows. The coefficient should be weighted based on the value of alpha1 in each cell, i.e. D = Dl*alpha1+Dg*(1-alpha1) where Dl and Dg are scalars. It should also include the addition of turbulent diffusion if alpha1 does not equal 0 or 1, i.e. not at an interface.

I included the following in the createFields.H header:

dimensionedScalar Dl
    (
        twoPhaseProperties.lookup("Dl")
    );

    dimensionedScalar Dg
    (
        twoPhaseProperties.lookup("Dg")
    );

I then created a volScalarField for D in createFields.H as :

Info<< "Reading field D\n" <<endl;
    volScalarField D
    (
      IOobject
      (
         "D",
         runTime.timeName(),
         mesh,
         IOobject::NO_READ,
         IOobject::AUTO_WRITE
      ),
       (Dl*alpha1+Dg*(scalar(1)-alpha1))
      );

Finally, I added the following to the CEqn.H file, which contains the transport equation:

volScalarField Dt = D+(turbulence->nut())/(scalar(0.7));

forAll(D, cellI)
   {
   if
   (
     alpha1[cellI]==scalar(1)
     ||
     alpha1[cellI]==scalar(0)
   )
      {
       D[cellI] = Dt[cellI];
      }
   }


    fvScalarMatrix CEqn
    (
         fvm::ddt(C)
         +fvm::div(phi,C)
         -fvm::laplacian(fvc::interpolate(D),C)
    );

    CEqn.solve();

Everything compiles and seems to run well. However, when I check the values of D, they aren't matching what I expect. In areas where alpha1=1 or 0, I am getting values of D that include turbulent diffusion. It's not clear to me why this is happening. Additionally, the D values shown in paraview look like they are not progressing in time. It appears they are adding to the existing D field each time step. However, the C values are clearly impacted by diffusion. What I want to do is change D each time step based on alpha1 values in each cell and write the results to the time directory to ensure I'm calculating correctly. Any advice on where I'm going wrong is appreciated.

Thanks,
Adam