Discussion Forum

February 9, 2012 1:54pm UTC

timeint (Time Integral) gives wrong answers for mass balance

Hello everyone!

I hope the COMSOL experts and enthusiasts here can help a fellow user out :)

The situation is that I have a tube (L=1m) with a porous medium in it and air flows in and out again and at the same time water is being evaporated from the porous solid phase to the gas phase. This 1D model uses "Species Transport In Porous Media" and "Darcy's Law". The first one calculates the concentrations, the second one gives the species mass and calculates the density and velocity. The 'drying' is modeled using a mass source of H2O that is dependent on temperature and the amount of water left (first order Arrhenius type equation). These I modeled myself and seem to be fine.

The problem is that I want to make a mass balance over the reactor. So, the time integral over all mass that flows in + the time integral over all mass that is released from a source term into the gas phase, is equal to the time integral over all mass that flows out at the exit. When the temperature is constant, this seems to be true with an 0.5% error, which is quite ok, but if the temperature changes over time (simulated using a step function) the difference in the time integral over what goes in and out is more than 50% off!! And so is the integral over the source flux, it is ridiculous. I may be doing something wrong here, but I have no idea what. Can anyone please help me?

I have included the file. It is not as complicated as it may seem. Just O2 and N2 (air) flow in without doing anything and there is 1 kg of water available for release. No water flows in, so the end point should 'see' 1kg of water come by over time. However, according to the timeint function the answer is 0.49 kg.. same for the integrated source flux.

Things I already tried:
1) Increasing the accuracy of the timeint function ( timeint(START,END,dl.rho*dl.u*A_flowthrough,1e-10) )
(1e-8 is standard), usually I use 1e-5 since it doesn't matter much and takes a lot longer to calculate).
2) fiddling around with the time stepping. Tried all kinds of settings, it converges well now, but still no good answer.
3) Changed Tolerances: Abstol is now 1e-3, Reltol is now 1e-8, quite okay I think.
4) Conservative vs Non-conservative formulation (compressible vs incompressible gas), does not matter much.
5) Changed amount of elements, currently 2000, just to be sure. But increasing # elements does not help.
6) Changed discretization order from Quadratic to Quartic and back again, no effect.

It must have something to do with temperature, because it did work when T was constant and now it doesn't since T is a function. Included in this post is the actual model (4.2a) and some graphs.
The first graph shows my own calculation of the amount of water left in the solid, this goes to 0, like it should.
The second one show the time integral over the exit and this is 0.49 kg while it should be 1 kg.
The third has the same problem.
The last one shows the H2O concentration in the entire reactor vs time. Most important is that is goes to zero, so apparently all mass flows out of the reactor.

I would be very happy with any useful input! Thanks in advance!

Reply  |  Reply with Quote  |  Send private message  |  Report Abuse

February 10, 2012 4:04pm UTC in response to Raymond Ruigbaard

Re: timeint (Time Integral) gives wrong answers for mass balance

Hi,

I always thought that the Darcy's law-Species transport in porous media combo was funny, here's why:

- If you solve for N species with Species transport in porous media, then you can readily calculate the pressure and density from ideal gas law (p=Ctot*R*T and rho=p*M/R/T), the fluid velocity from Darcy's law (calculate u=-K/mu*d(p,x) in a list of variables) --> you don't need to solve for the pressure

- the problem is overdefined when you solve for all N species and you use the Darcy's law interface, there is no need to solve for the pressure

Reply  |  Reply with Quote  |  Send private message  |  Report Abuse

February 10, 2012 6:22pm UTC in response to Nagi Elabbasi

Re: timeint (Time Integral) gives wrong answers for mass balance

sorry, I meant rho = sum(Ci*Mi) where Ci and Mi are species concentrations and molar masses.

Reply  |  Reply with Quote  |  Send private message  |  Report Abuse

February 12, 2012 3:03pm UTC in response to Raymond Ruigbaard

Re: timeint (Time Integral) gives wrong answers for mass balance

ni*Mi/V and Ci*Mi is not the ideal gas law, it's just a way to convert de values.

the sum of the species is the total mass of gas, so you already got the mass.

The Darcy's interface inserts the ideal gas law and the Darcy's law in the continuity equation, and the dependent variable of this eqn becomes the pressure. But you already calculated the pressure indirectly with the N Species Transport eqns : molar fraction % = partial pressure % = volume fraction % for an ideal gas and p = Ctot*R*T.


I didn't use 'timeint' when I tried the Darcy's law interface, but I had weird results : density calculated from the ideal gas didn't give the same result as rho=sum(Ci*Mi). I stopped using the Darcy's interface, only used the Species Transport and I calculated the density (mass concentration) and velocity in a list of variables, and I calculated the pressure p=Ctot*R*T with a coefficient form PDE because Comsol was complaining about a circular dependency if I put that equation in the same list of variables. This way the results are just fine and I validated my model by reproducing the results in a paper. It seems to be ok.

I think that the Darcy's law interface used with the Species Transport interface with N species is an overdefined problem and you get weird results. Maybe Comsol should clarify that.

Reply  |  Reply with Quote  |  Send private message  |  Report Abuse

February 13, 2012 2:29pm UTC in response to Raymond Ruigbaard

Re: timeint (Time Integral) gives wrong answers for mass balance

why can't you calculate the Darcy velocity? since you got the pressure from the species concentrations, you just have to calculate the velocity in a list of variables : u=-K/mu*d(p,x). The velocity will vary in the domain if the pressure varies.

If you absolutely want to use Darcy's interface, you have to set up the Species Transport interface with N-1 species... I made it work in the past, but the results were less good (numerical oscillations in the solution...).

Here's an example of how I set up my model. The attached model doesn't work at all, no boundary conditions adjusted, etc. It's just to show you where I calculate what, you'll understand my strategy to bypass the Darcy's interface.

enjoy!

Reply  |  Reply with Quote  |  Send private message  |  Report Abuse

February 15, 2012 11:30am UTC in response to Raymond Ruigbaard

Re: timeint (Time Integral) gives wrong answers for mass balance

Please?

Reply  |  Reply with Quote  |  Send private message  |  Report Abuse