Menu

Discussion Forum

New Discussion
Filter

Note: This discussion is about an older version of the COMSOL Multiphysics® software. The information provided may be out of date.

Discussion Closed This discussion was created more than 6 months ago and has been closed. To start a new discussion with a link back to this one, click here.

timeint (Time Integral) gives wrong answers for mass balance

Posted Feb 9, 2012, 8:54 a.m. EST Fluid & Heat, Chemical Reaction Engineering, Modeling Tools & Definitions, Parameters, Variables, & Functions, Results & Visualization Version 4.2a 11 Replies

Raymond Ruigbaard
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam


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!
11 Replies Last Post Feb 15, 2012, 10:53 a.m. EST
Raymond Ruigbaard
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 10, 2012, 8:01 a.m. EST
I said that with constant temperature I get 0.5% error but actually that is only when I exclude the Darcy's Law module, so only using Transport of species. If I include it again, I always get a 2% error (1.02 kg) for the time integration over the massflow at the exit boundary.

I have tried literally everything I can think of. Is this just the highest accuracy achievable? Is it because COMSOL uses FEM (finite element method) instead of FVM (finite volume method)? Am I doing something wrong?
Any thoughts here?
I said that with constant temperature I get 0.5% error but actually that is only when I exclude the Darcy's Law module, so only using Transport of species. If I include it again, I always get a 2% error (1.02 kg) for the time integration over the massflow at the exit boundary. I have tried literally everything I can think of. Is this just the highest accuracy achievable? Is it because COMSOL uses FEM (finite element method) instead of FVM (finite volume method)? Am I doing something wrong? Any thoughts here?
Francois Gregoire
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 10, 2012, 11:04 a.m. EST
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
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
Nagi Elabbasi Facebook Reality Labs
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 10, 2012, 12:17 p.m. EST
The two equations that you mentioned for ideal gas pressure and density are not independent. So, you cannot use them to separately evaluate both gas quantities.

Nagi Elabbasi
Veryst Engineering
The two equations that you mentioned for ideal gas pressure and density are not independent. So, you cannot use them to separately evaluate both gas quantities. Nagi Elabbasi Veryst Engineering
Francois Gregoire
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 10, 2012, 1:22 p.m. EST
sorry, I meant rho = sum(Ci*Mi) where Ci and Mi are species concentrations and molar masses.
sorry, I meant rho = sum(Ci*Mi) where Ci and Mi are species concentrations and molar masses.
Raymond Ruigbaard
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 11, 2012, 7:24 a.m. EST
I need both modules because I need to know the concentrations, the mass/density and the velocity. I can calculate the concentrations with the Species Transport module and from that I could calculate the density myself using Ideal Gas Law (rho_i=P_i*M_i/RT=n_i*M_i/V=C_i*M_i) but the velocity changes as well and for that I need Darcy. The problem is that Species Transport does not include mass, only species, so I mirrored the species in Darcy to make COMSOL give the species 'real' mass for mass balance purposes. For constant temperature this looks quite ok, as I said, but when temperature changes, the balance is way off.

Has anyone done mass balances on time dependent studies? (using timeint?) Or other, general ideas?
Thanks for any input.
I need both modules because I need to know the concentrations, the mass/density and the velocity. I can calculate the concentrations with the Species Transport module and from that I could calculate the density myself using Ideal Gas Law (rho_i=P_i*M_i/RT=n_i*M_i/V=C_i*M_i) but the velocity changes as well and for that I need Darcy. The problem is that Species Transport does not include mass, only species, so I mirrored the species in Darcy to make COMSOL give the species 'real' mass for mass balance purposes. For constant temperature this looks quite ok, as I said, but when temperature changes, the balance is way off. Has anyone done mass balances on time dependent studies? (using timeint?) Or other, general ideas? Thanks for any input.
Francois Gregoire
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 12, 2012, 10:03 a.m. EST
(rho_i=P_i*M_i/RT=n_i*M_i/V=C_i*M_i)
ni*Mi/V and Ci*Mi is not the ideal gas law, it's just a way to convert de values.

Species Transport does not include mass, only species, so I mirrored the species in Darcy to make COMSOL give the species 'real' mass for mass balance purposes
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.


Has anyone done mass balances on time dependent studies? (using timeint?)
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.
[QUOTE](rho_i=P_i*M_i/RT=n_i*M_i/V=C_i*M_i)[/QUOTE] ni*Mi/V and Ci*Mi is not the ideal gas law, it's just a way to convert de values. [QUOTE]Species Transport does not include mass, only species, so I mirrored the species in Darcy to make COMSOL give the species 'real' mass for mass balance purposes[/QUOTE] 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. [QUOTE]Has anyone done mass balances on time dependent studies? (using timeint?)[/QUOTE] 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.
Raymond Ruigbaard
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 13, 2012, 8:50 a.m. EST
the sum of the species is the total mass of gas, so you already got the mass.

Yes, it is indeed possible to calculate the mass yourself using the concentrations and the molar masses, so without using Darcy (and I did that with 0.5% error in the mass balance). But COMSOL does not understand that it is mass flowing through. I did a few simple tests and found out that if I use Species Transport AND Darcy together, without specifically assigning mass (in Darcy) the mass balance does not work since there is no mass flow.
But if I leave Darcy out, I can't calculate the changes in velocity (v=v_in everywhere then), at least not with the standard modules.


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.

That could be a problem, but the system does work as long as T is constant. If there was some conflict I would expect it not to work at all.

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 am curious how you managed to do that and would really like to see it. Would I ask too much if I asked if you could share that particular model? I won't be angry if you can't but I would be grateful if you can ;)


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.

Interesting idea, I think I'll ask them. Thanks for the input.

[quote]the sum of the species is the total mass of gas, so you already got the mass.[/quote] Yes, it is indeed possible to calculate the mass yourself using the concentrations and the molar masses, so without using Darcy (and I did that with 0.5% error in the mass balance). But COMSOL does not understand that it is mass flowing through. I did a few simple tests and found out that if I use Species Transport AND Darcy together, without specifically assigning mass (in Darcy) the mass balance does not work since there is no mass flow. But if I leave Darcy out, I can't calculate the changes in velocity (v=v_in everywhere then), at least not with the standard modules. [quote] 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. [/quote] That could be a problem, but the system does work as long as T is constant. If there was some conflict I would expect it not to work at all. [QUOTE] 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. [/quote] I am curious how you managed to do that and would really like to see it. Would I ask too much if I asked if you could share that particular model? I won't be angry if you can't but I would be grateful if you can ;) [quote] 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. [/QUOTE] Interesting idea, I think I'll ask them. Thanks for the input.
Francois Gregoire
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 13, 2012, 9:29 a.m. EST
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!

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!

Attachments:

Raymond Ruigbaard
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 14, 2012, 6:25 a.m. EST
Dear Francois,

I implemented your method but unfortunately it is now not solving at all. There seem to be convergence issues. Maybe you would be so kind to take a quick look to see where my model is different from your (working) model?
Maybe I overlooked some BC or something. Also, I know that Darcy works best with relative pressure since the differences are very small. This method works with the absolute pressure right? Could there be problems there?

Thanks for the help and hope to hear from you soon :)

Ray.

(Current model included)
Dear Francois, I implemented your method but unfortunately it is now not solving at all. There seem to be convergence issues. Maybe you would be so kind to take a quick look to see where my model is different from your (working) model? Maybe I overlooked some BC or something. Also, I know that Darcy works best with relative pressure since the differences are very small. This method works with the absolute pressure right? Could there be problems there? Thanks for the help and hope to hear from you soon :) Ray. (Current model included)
Raymond Ruigbaard
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 15, 2012, 6:30 a.m. EST
Please?
Please?
Francois Gregoire
Send Private Message Flag post as spam

Please login with a confirmed email address before reporting spam

Posted: 1 decade ago Feb 15, 2012, 10:53 a.m. EST
Hi Ray,

Sorry I don't have the time to analyse your model in details and make it converge. But I may have some remarks:

1. in the inflow equations ( ex.: (1*(P+p)*(1*V_fraction_N2))/(R*T_S) ), why do you multiply by 1 at the numerator? and why p = 0 in the parameters list? this expression could be (P*V_fraction_N2)/(R*T_S)

2. you should apply a pressure dirichlet boundary condition in the PDE node for ideal gas law. The pressure applied at the inflow should the same as the one underlied in the inflow conditions in Species transport

3. the sum of the diffusive fluxes of the 4 species should be equal to 0 in order to theoretically find back the continuity equation when you sum all species equation (read theory about multicomponent species transport). So you should modify the diffusive flux equation of 1 species in order to respect that. ex.: dfluxx_N2 = -dfluxx_O2 - dfluxx_H2O - dfluxx_CO2

4. use a coarser mesh to make things quicker, no need to have a super mesh at first.

Hope it works
Hi Ray, Sorry I don't have the time to analyse your model in details and make it converge. But I may have some remarks: 1. in the inflow equations ( ex.: (1*(P+p)*(1*V_fraction_N2))/(R*T_S) ), why do you multiply by 1 at the numerator? and why p = 0 in the parameters list? this expression could be (P*V_fraction_N2)/(R*T_S) 2. you should apply a pressure dirichlet boundary condition in the PDE node for ideal gas law. The pressure applied at the inflow should the same as the one underlied in the inflow conditions in Species transport 3. the sum of the diffusive fluxes of the 4 species should be equal to 0 in order to theoretically find back the continuity equation when you sum all species equation (read theory about multicomponent species transport). So you should modify the diffusive flux equation of 1 species in order to respect that. ex.: dfluxx_N2 = -dfluxx_O2 - dfluxx_H2O - dfluxx_CO2 4. use a coarser mesh to make things quicker, no need to have a super mesh at first. Hope it works

Note that while COMSOL employees may participate in the discussion forum, COMSOL® software users who are on-subscription should submit their questions via the Support Center for a more comprehensive response from the Technical Support team.

Suggested Content