Finite difference transport model

[Petrusvspreng]

I am running a TRANSPORT model as per example 13C, using MIX equations to simulate diffusion in spherical coordinates. I use 20 vertical increments and 5 immobile cells per vertical increment.

I have 1.15 mole/kg_w of a kinetic compound (CuO) in each of the immobile cells (zero in mobile cells). Each cell (1-121) with 1 kg water.

A sulphuric acid solution is percolated via solution 0 which dissolves the CuO to yield Cu+2 into solution.

I keep to the relation between irrigation rate, shift duration and vertical increment to ensure the solution contained in one vertical increment is moved through the mobile cells per shift.

The resulting profile of Cu+2 concentration in the drainage solution (20), and the extent of Cu extraction based on Cu+2 in drainage solution check virtually perfectly with other well established software against which I calibrate it. (I count the moles of Cu+2 passing via solution 20 with each shift using the TotMole("Cu") variable, and its molality using Tot("Cu")  ).

But then I thought to go one step further and check the mass balance between solid and solution phases.  I count the remaining moles of CuO in solid phase using KIN("CuO") and as another check use KIN_DELTA("CuO") to check the amount of CuO having reacted. To my surprise I find that more moles of Cu+2 are passing via solution 20 than moles of CuO that have disappeared from the solid phase. In terms of actual numbers, I find that of the 100x1.15=115 moles of CuO I started with, 80 moles passed via solution 20 in 150 days, yet only 42 moles of CuO have disappeared from the solid phase. There is an error by factor of 2, it would make sense if 80 moles passed as drainage while 84 moles of CuO has disappeared from the solids, the difference of 4 moles would be accounted for by Cu+2 in the solution phase that has not yet passed the drainage point.

I'll probably kick myself when I realise my mistake, but for now it has me stumped.

- Petrus van Staden.

[dlparkhurst]

Maybe I need to see the input file for this one.

[Petrusvspreng]

OK find it attached as a PQI.

I normally run it from Excel, relying on the COM, which directs PUNCH'es to my spreadsheet from where I plot and calculate. The graph produced during run time is not very useful other than indicating progress towards completion.

For your purposes I have set the TRANSPORT block to run only 300 shifts because it is not the fastest.

[dlparkhurst]

Are you sure your mixing fractions are correct?  My guess is that you are losing/gaining Cu because of the mixing fractions.

Try simplifying your system to one cell with one stagnant cell or column to track down the problem.

[Petrusvspreng]

Well my calcs replicate the Example 13C mixing factors exactly for the conditions in that example, after setting up my own calcs to determine the Vj, Aij, h, fbc etc and from there the mixing factors.

My case only deviates by considering different values of radius and mobile/immobile solution fractions, which therefore yield mixing factors that differ from those of Example 13C but using the same calculation logic.

Any case I will revert to the simpler first-order exchange version and see how that works, thanks.

[dlparkhurst]

Does example 13C lose mass by your calculations?

[Petrusvspreng]

Thing is Example 13C does not use chemical reaction with a solid phase. I know that it does maintain a mass balance between soluble species in the in- and outlet (I corresponded with you about that previously when I had an apparent +/-15% mass balance error between in and out, but that turned out to be as a result of my simulation involving a pulse test with very concentrated inlet pulse, I was calculating the balance based on shift volumes, where I was supposed to work with kg water per shift, and with such concentrated inlet solution there is a big difference between molarity and molality).

What I have not done is to check a dynamic balance for Example 13C between all of what had gone in, out, and retained at any given time. But even if that checks out it does not speak directly to my example.

For my example I need to see a balance between (a) moles of kinetic (solid phase) reactant in all immobile cells at time=0 and (b) moles of reactant remaining in immobile cells PLUS what has passed out the bottom at time = t (plus a bit that will have gone into solution but still remains in the column which is a relatively small amount). That does not balance and there is not an example like that which I know of. And of course the weirdest thing is that more seems to have passed out the bottom than has dissolved from the solid phase. So I am not actually losing copper, I am gaining it.

[dlparkhurst]

13C has EXCHANGE, which is a solid phase, so a mass balance on an exchangeable cation should be analogous.

Seems like you may need to take into account the different volumes of the stagnant cells. By summing the KINETIC reactants in each cell equally, you may be assuming that the volumes of each cell are the same, whereas the stagnant cells are smaller than the mobile cell, and each stagnant cell volume varies as it approaches the center of the sphere.

Clearly, I am not very good with the mixing fraction approach. Tony is the person you probably need to ask.

[Petrusvspreng]

I have been wondering about that, but then discarded the thought because the programme does not 'know' that cell number # represents a cell towards the centre or outside. From the theoretical point of view the mixing factors are supposed to represent the rate at which leaching will take place from a sphere, but with equi-sized immobile cells, which is indeed the way in which Ex13C was set up.

But I will check it out nevertheless. The outermost cell already accounts for almost 50% of the sphere volume and of course sees most of the leaching taking place. Could  just be where my factor of 2 is hiding.

I seem to be receiving a disproportionate share of your efforts but its good to be able to throw some thoughts around, thanks.

[Petrusvspreng]

Hi,

Just as a follow-up to the discussion on this topic. I have now done the first-order exchange model as well and the principle I have figured out is:
a. Sure the SOLUTION blocks define all solutions on the basis of 1kg of water (at least by default) and the MIX'es of the finite difference diffusion system (Example 13C) are set up using equi-sized solutions to emulate the correct concentration-vs-time response,
however
b. When it comes to the mass balance, one needs to consider the characteristic of the actual physical system being modelled. That involves that one selects a good reference solution and in the mass balance calculations do the necessary proportioning of the other solutions. For example I can select as a basis the mass of mobile solution that moves per shift. From what I know about the physical system I need to proportion the amount of solution considered to be present in each immobile cell, the same with the solid phase(s) present. Parameters like KIN_DELTA indicates what happened per 1kg of solution basis, but mass balance calc's need to proportion that for the actual kg amount of that solution being present based on the calculation basis selected.

 

[dlparkhurst]

I hope your balances are adding up now. Sorry I could not be more help.

[Petrusvspreng]

Things do balance now and it was your suggestions that put me on the right track, thanks.