[Rafael Jurado]

Hello everyone,

I would like to modify the IcoFoam solver to solve two MassTransport equation which are coupled like this:

1.- ∂C1/∂t + ∂Ui·C1/∂xj = D11 · ∂^2 C1/∂^2xj + D12 · ∂^2 C2/∂^2xj
2.- ∂C2/∂t + ∂Ui·C2/∂xj = D22 · ∂^2 C2/∂^2xj + D21 · ∂^2 C1/∂^2xj

I would like to do this time discretization:

(C1^(n+1) - C1^(n))/dt = 0.5*[DiffC1 + DiffC2]^(n+1)
+ 0.5*[DiffC1 + DiffC2]^(n)
+ 1.5*ConvC1^(n) - 0.5*ConvC1^(n-1)


(C2^(n+1) - C2^(n))/dt = 0.5*[DiffC2 + DiffC1]^(n+1)
+ 0.5*[DiffC2 + DiffC1]^(n)
+ 1.5*ConvC2^(n) - 0.5*ConvC2^(n-1)

So, I have tried this in OpenFOAM:

1.-
solve
(
fvm::ddt(C1)
+ 0.5*(fvm::div(phi, C1) + fvc::div(phi, C1))
- 0.5*(fvm::laplacian(D11, C1) + fvc::laplacian(D11, C1))
- 0.5*(fvm::laplacian(D12, C2) + fvc::laplacian(D12, C2))
);

2.-
solve
(
fvm::ddt(C2)
+ 0.5*(fvm::div(phi, C2) + fvc::div(phi, C2))
- 0.5*(fvm::laplacian(D22, C2) + fvc::laplacian(D22, C2))
- 0.5*(fvm::laplacian(D21, C1) + fvc::laplacian(D21, C1))
);

and I have putted this into a loop to iterate until C1 and C2 is converged.

The problem comes from the term [fvm::laplacian(D12, C2)] in the first equation and [fvm::laplacian(D21, C1)] in the second one. I think that this is because, you can not use an implicit term of a variable while solving the other one. Let say the term [fvm::laplacian(D12, C2)] (i.e. DiffC2^(n+1)) in the first equation.

So I have two question, the first one and the most important is:
1.- How can I managed to introduce the other Concentration explicitly in each equation in OpenFOAM?

And the second one is;

2.- How can I take the concentration variable of two previous timesteps for using an Adam-Bashforth time discretization in the convection terms, i.e. [-0.5*ConvCi^(n-1)]?

Thank you in advance...

Kind regards,

[Rafael Jurado]

Anyone can help me please?

Thank's in advance

[chriss85]

If you use fvc instead of fvm you can calculate the laplacian explicitly. The fvm operators all concern the variable for which you solve the equation, so using the fvc operators uses the fields from the last time that variable was solved for (or modified otherwise).

[Rafael Jurado]

Thanks for your reply chriss85,

But the problem is that I want to if it is possible to use the "current" timestep value of the component, let say C2 in C1 equation, in order to use a second order scheme such as Crank-Nicolson. So for doing it I need the (n+1) time of C2... I don't know if I have explain it porperly...

Thanks you again.

[chriss85]

I think you can use field.oldTime() and possibly field.oldTime().oldTime() to access the old fields from previous time steps, but I'm not too sure about this. If this doesn't work you'll have to copy your old fields somewhere appropriate. If you need a different time derivative discretization you'll have to write a class for this. In this case you should look up the implementation of the existing time derivative operator implementations for reference.

[Rafael Jurado]

I have try this already, but it doesn't work with the field.oldTime(). So, in simple words:

- Is there any way to update the value of C1 and C2 within the same timestep?

If I could do this, I think that it would be enough with put the two equations into a loop...

Thanks for your quick reply.