[gergelyacs]

Hello Foamers,

I'm developing a solver for electrochemical flow simulations, and I have to solve an equation, which has the general form (similar to the incompressible continuity equation):
div(flux)=0, which I'll to solve as

fvScalarMatrix UeEqn
(
fvc::div(Flux)
);

where Flux is a volVectorField. However, when I try to compile my solver, I always get error messages, I've tried also some other ways to code it, but without success so far.
I would be glad if you could help me with this problem.

Regards,
Gergely

[chriss85]

You will likely need a reformulation such as div(grad(phi))=0 and then calculate flux=grad(phi) or so.

An example for this is magnetostatics with a constant current, phi is electrostatic potential, and what you want to solve is div(j)=0 with j being current density (j=-sigma*grad(phi)). What you end up solving is laplacian(sigma, phi)=0, i.e. div(sigma*grad(phi))=0. Also keep in mind that div(flux)=0 doesn't have a unique solution.

[gergelyacs]

Thank you, chriss85, indeed expanding the equations did help, my solver could be compiled, however, it lead me only to "incompatible fields", thus I expanded my equation even deeper.
My eqaution is
div( (D1-D2)*grad(C)+(D1+D2)*F/R/T*C*grad(Ue) ).
Where D1,D2,F,R & T are dimensioned scalar constants, C and Ue are my variables, both are scalar fields.
After expansion I've found:
laplacian( (D1-D2)*grad(C))+ (D1+D2)*F/R/T* ( grad(C)*grad(Ue)+C*laplacian(Ue) ).
Translated into OpenFOAM:
fvm::laplacian(D1-D2,C)+ (D1+D2)*F/R/T* (fvc::grad(C,Ue)+C*fvc::laplacian(Ue))

In my understading it should be correct, however the "fvc::grad(C,Ue)" terms with the following error message:
"ElectricEqn.H:8:16: note: ‘Foam::volScalarField {aka Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>}’ is not derived from ‘const Foam::tmp<Foam::GeometricField<Type, Foam::fvsPatchField, Foam::surfaceMesh> >’
fvc::grad(C,Ue)"

I'm really not sure what does that mean, and how could it be solved. I would be glad for any help.

[chriss85]

I'm not quite sure where the problem comes from without seeing the type definitions of the variables.

Some ideas:
Are you trying to solve for both C and Ue? If so, you will need a second equation for Ue. Also keep in mind that by explicitly using C in the equation the solver will use the field from the last iteration there.

A tmp<> is basically a pointer to the object to save performance by avoiding copies at assignments. You can dereference it with (), which is sometimes needed.

[gergelyacs]

Thank you for the suggestions.
Yes, I have two equations, one is the previously discussed transport equation, the other is a "classic" Nernst-Planck equation, where the time-derivative of the concentration equals the gradient of the concentration and the potential field. With that I don't have any problems.
It is actually very similar to the other equation, the main difference is that the time-derivative is missing due to some analytical reasons.
Lately I've made a formulation which will enable my code to run, where I've separated my equation in the following way:

FluxUe=(D1+D2)*F/R/T*C*fvc::grad(Ue);
fvScalarMatrix UeEqn
(
fvm::laplacian(D1-D2,C)+fvc::div(FluxUe)
);
However, this fails to run (Floating point exception) wherever my potential field is not uniform, thus I fear this formulation is essentially incorrect.
I'm attaching my code with a very simple 1D test-case, if that would help. In this state may it doesn't makes sense, why would I use OpenFoam, but the next step will be to couple this system with the incompressible Navier-Stokes equations.

[chriss85]

In which of these two lines is it failing?
Have you checked for division by zero in the first line?

[gergelyacs]

Actually none of these two lines are failing but when I call the solve command with "UeEqn.solve();", only in such a case when the first line creates a non-zero field.
The divisions in the first line are only constant positive scalars, and the whole FluxUe field had only finite values, in other words there wasn't any NaN or Inf values.