[pbe_cfd]

Hi,

I need to transport n (>>1) active scalar within Oper. Split. Approach. I have a code to solve the reaction step so, what I do is first solving the reaction step for cell based values then solving the transport step. A pseudo code is given here

1. Step: Solve reaction for each cell

for (i=0; i<C[0].size(); i++) // C := PtrList<volScalarField>
// Assign values of all the scalars at the i'th cell to Ci
for ( j=0; j<n; j++)
{
Ci[j]=C[j].internalField()[i]; Ci: is an array of length n, for the values of n Scalar at the i'th cell
}
solveReaction(Ci) // Ci: is an in/output and updated in the routine solveReaction
for ( j=0; j<n; j++)
{
C[j].internalField()[i]=Ci[j]=; Updating the values of C at the i'th cell with the updated values of Ci
}

end

2. Step: Transport each scalar field, C[i]
for (i=0; i<n; i++)
{
solve
(
fvm::ddt(C[i])
+ fvm::div(phi, C[i])
- fvm::laplacian(nuT, C[i]),
mesh.solver("class")
);
}

I implemented this pseudo code by extending scalarTranportFoam. The problem is, although each step works fine independently:
1) Cancel 2. step and update the scalar fields by reaction. Result: OK
2) Cancel 1. step and solve the passive transport problem for n Scalar. Result: OK

When I calculate them together, the results are as if I had not performed the reaction step. What can be the problem??? And, how can it be fixed ???

Thanks in advance,
Evren

PS. I can post or can send the piece of code via e-mail to whom are interested in.

PPS. Boundary values are another problem which is hacked by using symmetryPlane BC for the moment.