Dear Friends, I am trying to solve a electromigation problem. I have only diffusion active and convection inactive in terms. I am trying to add the migration term as source. The migration term consists of gradient of (u*c*(delPhi/delx)).Means I am calculating the migration fluxes at each face (east and west), subracting west value from east to get source or sink due to migration in that cell. Since c and phi (volts) are scalar, I am calculaing the values of c and delphi/delx at wall faces by taking average cell vales of c and delphi/delx of neighboring cells. Now my problems are (i) Solution is not converging and the values of c calculation by this route are oscillating. Is there any command in phoenics to calculate the mass fluxes at walls? (ii) at boundaries, ie 1 and nx grid, how do I determine the c and delphi/delx values at west and east faces respectively. As outside boundaries I do not have phi exising.

Thanks in advance

Regards, Sanjay

In continuation from above problem, in fact if can calculate the scalars at walls (not leading to inconsistancy) such as c and delphi/delx from the cell center values, please let me know. Is the arthematic averaging I am doing to find the scalars at wall faces using cell centre values advisable? Regards, Sanjay

What is the dependent variable of the conservation equation containing the migration term? What is c and what is u?

Dear Mike My Equation consists of del(c)/del(t) = -D(del2(c)/del(x)2)+del(u*c*del(V)/del(x)) where c=concentration

D=diffusive coefficient

the last term on RHS is migration term where u=mobility, V=volts. Using PHOENICS, I am solving this equation, Diffusion goes straight forward as -prandl, I do not want to solve the migration term (which is actually convective ion motion) by convection as there will be problems with continuity. The migration is for ions, not fluid. So I want to solve the migration term by source/sink in each cell. So my intention is to find the u and del(V)/del(x) at walls from the cell centre values available for them and then subtract the flux of east face from west face and introduce the same as fixflu in patch/coval [coval(acsoil,c3,FIXFLU,grnd)] (source solving by group 13). My problems here are that I am doubtful whether I can subtract the fluxes of boundaries and introduce the same as source is advisable or not. As source and sink occurs (since del (V)/del (x) varies in spacial coordinate) there is possibility of the difference of fluxes calculated will have either positive or negative sign. Does it mean that the flux is changing the direction wich actually means source or sink depending on the sign? any how my solution is oscillating. Can someone suggest how to change the coefficients of algebric equation after discritisaton in PHOENICS (ie. ap, ae, aw)?

Thanks, now it is clear to me what you are trying to do. This is what I suggest: interpret -u*dphi/dxi as a "migration" velocity vi so that the source term in the phi equation is -d/dxi(vi*c). The source term can then be represented as sum(migration flux inflow)-sum(migration flux outflow), where the summation is over all faces of the cell, but you upwind interpolation must be used for each contribution. In other words, the source term should be discretised like a convection term, rather than as a diffusion term, otherwise numerical instabilities will result.

Thank you dear Mike! Let me try and I will let you know the result! Regards, Sanjay

Dear Mike, I understood your interpretation, but please let me know if I have to write the code for migration using upwind scheme in Group13 or is there a way out in PHOENICS? I hope with convection off and writing fortran code for upwind scheme for C using velocity vi at walls will subside the oscillations. In that case can I calculate the velocities at walls by averaging cell values of -u*dphi/dxi and choose c as per upwind scheme? Thank you for all the efforts you are taking. Regards, Sanjay

Dear Mike, I used as per your suggestion I tried to implement in two different ways (1)My discritisation is boiling down to uCw(dphi/dx)e-uCp(dphi/dx)w using upwinding scheme, where Cw is concentration in west cell, d(phi/dx)e is potential derivative at east wall and u is the mobility etc. I calculated from IX=1,nx (since I am solving 1-Dimensional problem only) the fluxes and introduced as net of above equation in patch/coval through Group 13 as

patch(acsoil,CELL,1,nx,1,1,1,1,1,lstep)

coval(acsoil,c3,fixflu,grnd2)

where grnd2 consists of the net value of equation above. Here the net value is either positive or negative which is expected as source and sink occurs in different regions. So the solution is oscillating. Then I changed the patch and covals to introduce the Co and Val as suggested by you.... (2) The patch and coval are

patch(acsoil,CELL,1,nx,1,1,1,1,1,lstep)

coval(acsoil,c3,grnd2,grnd2)

and in Group 13, the corresponding Co is -u(dphi/dx) which is always negative and val is uCw(dphi/dx)e.

There is drastic upward change in the all values now with this scheme. There are no oscillations, but as the values are shooting up every node. Is this the way to implement convective sources? Am I correctly putting the co and val values? Please help me. Thanks in advance! Regards, Sanjay

The source term should be implemented in group 13, and the "conventional" convection terms need only to be switched off if they do not appear in the conservation equation for the scalar variable C. The migration velocity should be computed at the cell faces, which requires no averaging for dphi/dx, but u will be averaged from the cell nodes to the cell face. The migration velocities will be upwinded at the cell faces, which then determines the neighbour value of c to be convected across the cell face by the migration velocity.

The precise form of the source term depends on the sign of the migration velocity at each face. I don't think I actually suggested how to implement the source term through CO and VAL. I haven't the time to think this through but I expect you to have two terms, one for each face, and only one of these can be linearised and the other cannot. Therefore, linearisation is possible only for that term which involves outflow of the nodal value of the scalar, and the term involving inflow of the neighbour scalar value must be introduced explicitly as a source term. I hope I have got the signs right, but you can work out the detail. Certainly, CO cannot be allowed to become negative.

Dear Mike, Thank you. I will try. regards, Sanjay