CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Running, Solving & CFD (
-   -   Second equation uses wrong state variables boundary condition (

Fabian_W January 22, 2013 08:08

Second equation uses wrong state variables boundary condition

this is my first post, so I introduce myself quickly. I am Fabian, a Postdoc in Germany. I do topology optimization and now work on an OF based project. Someone else wrote most of the solver and other stuff, my task is to change the code in such a way, that a second problem is solved afterwards. I do it in the same code.

I have good C++ knowledge but am totally new to OpenFOAM.

The forward problem is solved for h with dynamic boundary conditions defined based on the time in 0/h as inlet code. The forward problem is nonlinear with implicit time stepping.

I need to solve a so called adjoint problem which is linear, needs to be solved with explicit time stepping and goes backwards in time. I do the backwards in time by going forward in time (resetting my time object) but have for the solution variable the boundary conditions with other times in the dynamic 0/lambda code.

I have

fvScalarMatrix lambda_eqn(
phi * C * fvm::ddt(lambda) == -fvc::laplacian(k, lambda)

and hope, this solves my problem with an explicit time stepping.

in fvSolution I have

solver PCG;
preconditioner DIC;
// tolerance 1e-06;
tolerance 1e-14;
relTol 0.1;

solver PCG;
preconditioner DIC;
// tolerance 1e-06;
tolerance 1e-14;
relTol 0.1;

relTol 0.0;

In fvSchemes I have

default none;
laplacian(k,lambda) Gauss linear uncorrected;
laplacian(k,h) Gauss linear uncorrected;

For the forward problem the inlet code of 0/h is executed (I do a command line output) and then it is written "DICPCG: Solving for h, Initial residual ..."

For the adjoint problem again the code from 0/h is executed and then it is written " diagonal: Solving for lambda, Initial residual = 0, Final residual = 0, No Iterations 0"

My problem is, that the inlet code from h is executed but not from lambda ?!

Any help or hint is very much appreciated!

Thank you very much, Fabian

chegdan January 22, 2013 19:18


there are a number of things that could be the issue. You could...

* Have a zero values for your constants that zero out your field
* Ill-posed boundary conditions (e.g. all are set to zeroGradient) that essentially do nothing
* Mean to use fvm::laplacian(k,lambda) instead of fvc::laplacian(
* Have an issue somewhere else in your algorithm that is causing your problem.

is there more code that you have changed that we can see and track down the problem?

Fabian_W January 23, 2013 05:00


thank you very much for reply. Please note, that I am really new to OF and have much
more guessing about how it works than knowledge.

From you points I can see that they are issues, but I do not see, how they could cause my problem, that the boundary conditions from the state variable of the first problem are used for the second problem with another state variable.


Originally Posted by chegdan (Post 403439)
* Mean to use fvm::laplacian(k,lambda) instead of fvc::laplacian(

As far as I understood this is the notation for an explicit in time scheme. As there is no explicit ddt scheme I could find, I hope that OF uses explicit from my notation.

Here is some of the code

Foam::Time time(Foam::Time::controlDictName, args);
Foam::fvMesh mesh(Foam::IOobject(Foam::fvMesh::defaultRegion, time.timeName(), time, Foam::IOobject::MUST_READ));

// h is state variable for first problem
volScalarField h(IOobject("h",time.timeName(),mesh,IOobject::MUST _READ,IOobject::AUTO_WRITE),mesh);
// state variable for second problem (adjoint problem) adjoint solution - counterpart of h. Reversed time scheme in 0/lambda
volScalarField lambda(IOobject("lambda",time.timeName(),mesh,IOob ject::MUST_READ,IOobject::AUTO_WRITE),mesh);

volVectorField rho(IOobject("rho",time.timeName(),mesh,IOobject:: MUST_READ_IF_MODIFIED,IOobject::AUTO_WRITE),mesh);
volScalarField phi(IOobject("phi",time.timeName(),mesh,IOobject:: READ_IF_PRESENT,IOobject::NO_WRITE),mesh);
volScalarField s(IOobject("s",time.timeName(),mesh,IOobject::MUST _READ,IOobject::AUTO_WRITE),mesh);
volScalarField k(IOobject("k",time.timeName(),mesh,IOobject::MUST _READ,IOobject::AUTO_WRITE),mesh);
volScalarField C(IOobject("C",time.timeName(),mesh,IOobject::MUST _READ,IOobject::AUTO_WRITE),mesh);


fvScalarMatrix hEqn
phi*C*(fvm::ddt(h) - fvc::ddt(h)) == fvm::laplacian(k,h) - phi*fvc::ddt(s)

// now adjoint problem
time.setTime(0.0, 0); // reset to start again from 0

// history saves h for every timestep in the forward problem
for (int i = history.size() - 1; i >= 0; i--)
const TimeStep& step = history[i];

// set C and K
fvScalarMatrix lambda_eqn(
-1.0 * phi * C * fvm::ddt(lambda) == -fvc::laplacian(k, lamba)

Now the problem is, that during solving the first and second problem, the dynamic boundary conditions from the inlet code 0/h are used but not from 0/lambda. But 0/lambda is compiled (easy to check with a syntax error)!

Why is OF touching h in the second problem?

Thanks for any help!

All times are GMT -4. The time now is 01:10.