Re: Solving advection, segregation and diffusion equation

[James Pringle]

Dear all -- in a bit I am going to post a little idealized
advection/diffusion problem to give an idea of how some of these errors in
steady solutions arise, and how iterating in time helps. It will be much
more basic then this work, but you might find it helpful. It was very
helpful for me in figuring out my more complex system to reproduce the
error in a simple way. It also reveals some surprising issues in the fipy
code -- for example, the error in a steady solution can depend on the
"initial" condition you feed the solver.

It will be under the title "sources errors in advection/diffusion problems,
and one solution." I hope to post it tonight.

Jamie

On Thu, Sep 15, 2016 at 1:27 PM, Zhekai Deng <
zhekaideng2...@u.northwestern.edu> wrote:

> Thanks for the reply. Those are very useful suggestions. I have some new
> findings which may help shedding some light on this problem.
>
> I have added TransientTerm and tried to see how the solution is developed
> over time. I also remove the datamin and datamax input on the viewer, so
> that I can see the full range of the solution variable data.
>
> When I see the temporal development of the solution (attached in the
> email), It looks stable. However, I noticed that the right outlet boundary
> has extreme high concentration and did not outflow the flow domain. I think
> It has due to the improper setup of my right outlet condition. Since the
> materials keep building up on the right outlet boundary, it appears to
> affect the solution in the actual flow domain(also attached in email) and
> drives the solution unstable. This could be the main reason that solving
> steady state equation directly never converges.
>
> I then begin to look into how to properly set up "open outlet" boundary,
> the closest I could find is following:
>
> https://www.mail-archive.com/fipy@nist.gov/msg02797.html
> 
> and in the section of "Applying fixed flux boundary conditions" of
> official website
> http://www.ctcms.nist.gov/fipy/documentation/USAGE.html
> 
>
> I tried to adapt the official method listed on the official website, but
> no success is achieved.  Thus, I wonder is it possible to point out some
> examples that used this method to define fixed flux outlet boundary
> condition ? You may find my current progress on this part in the code I
> attached.
>
> Another thing I noticed is that,  there is a reported issue on robin
> boundary condition:
>
> https://github.com/usnistgov/fipy/issues/426
> 
>
> Would it be possible to explain more about this ? I don't think the up &
> bottom boundary condition in my case causes the non convergent problem,
> because from the Transient equation, the segregation of the material starts
> to develop at the top&bottom boundary, which is consistent with my
> expectation.
>
> I have attached my current progress of the code, and the temporal
> development of the solution.
>
> Best,
>
> Zhekai
>
>
>
>
>
>
> On Thu, Sep 15, 2016 at 9:33 AM, Daniel Wheeler  > wrote:
>
>> On Wed, Sep 14, 2016 at 7:52 PM, Guyer, Jonathan E. Dr. (Fed)
>>  wrote:
>> >
>> > I don't know offhand. With a Péclet number of 100, your problem is
>> almost completely hyperbolic, which FiPy (and cell-centered Finite Volume)
>> isn't very good at. Daniel knows more about this and may have some
>> suggestions.
>> >
>> > You might consider adding a TransientTerm for artificial relaxation and
>> trying the VanLeerConvectionTerm.
>>
>> Definitely use a TransientTerm and see the time step you can get away
>> with. Don't sweep, just use time steps. Using the VanLeerConvetionTerm
>> isn't necessary as the accuracy is not important in a steady state
>> problem.
>>
>> --
>> Daniel Wheeler
>>
>> ___
>> fipy mailing list
>> fipy@nist.gov
>> http://www.ctcms.nist.gov/fipy
>> 
>>   [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy
>> 

Re: Solving advection, segregation and diffusion equation

[Zhekai Deng]

Thanks for the reply. Those are very useful suggestions. I have some new
findings which may help shedding some light on this problem.

I have added TransientTerm and tried to see how the solution is developed
over time. I also remove the datamin and datamax input on the viewer, so
that I can see the full range of the solution variable data.

When I see the temporal development of the solution (attached in the
email), It looks stable. However, I noticed that the right outlet boundary
has extreme high concentration and did not outflow the flow domain. I think
It has due to the improper setup of my right outlet condition. Since the
materials keep building up on the right outlet boundary, it appears to
affect the solution in the actual flow domain(also attached in email) and
drives the solution unstable. This could be the main reason that solving
steady state equation directly never converges.

I then begin to look into how to properly set up "open outlet" boundary,
the closest I could find is following:

https://www.mail-archive.com/fipy@nist.gov/msg02797.html
and in the section of "Applying fixed flux boundary conditions" of official
website
http://www.ctcms.nist.gov/fipy/documentation/USAGE.html

I tried to adapt the official method listed on the official website, but no
success is achieved.  Thus, I wonder is it possible to point out some
examples that used this method to define fixed flux outlet boundary
condition ? You may find my current progress on this part in the code I
attached.

Another thing I noticed is that,  there is a reported issue on robin
boundary condition:

https://github.com/usnistgov/fipy/issues/426

Would it be possible to explain more about this ? I don't think the up &
bottom boundary condition in my case causes the non convergent problem,
because from the Transient equation, the segregation of the material starts
to develop at the top&bottom boundary, which is consistent with my
expectation.

I have attached my current progress of the code, and the temporal
development of the solution.

Best,

Zhekai






On Thu, Sep 15, 2016 at 9:33 AM, Daniel Wheeler 
wrote:

> On Wed, Sep 14, 2016 at 7:52 PM, Guyer, Jonathan E. Dr. (Fed)
>  wrote:
> >
> > I don't know offhand. With a Péclet number of 100, your problem is
> almost completely hyperbolic, which FiPy (and cell-centered Finite Volume)
> isn't very good at. Daniel knows more about this and may have some
> suggestions.
> >
> > You might consider adding a TransientTerm for artificial relaxation and
> trying the VanLeerConvectionTerm.
>
> Definitely use a TransientTerm and see the time step you can get away
> with. Don't sweep, just use time steps. Using the VanLeerConvetionTerm
> isn't necessary as the accuracy is not important in a steady state
> problem.
>
> --
> Daniel Wheeler
>
> ___
> fipy mailing list
> fipy@nist.gov
> http://www.ctcms.nist.gov/fipy
>   [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
>
from fipy import *

# mesh set up
nx = 100.
nz = 100.
dx = 1/nx
dz = 1/nz

mesh = Grid2D(dx=dx, dy=dz, nx=nx, ny=nz)
phi = CellVariable(name = "solution variable",mesh = mesh, value = 0.5)


# the follwing set up the velocity and shear rate in vector form.

# for those direction that I don't need it, I set its value as zero.

velocityX = mesh.faceCenters[1] # linear velocity profile
velocityY = FaceVariable(mesh=mesh, value = 0) # zero velocity in y direction
velocityVector = FaceVariable(mesh=mesh, rank=1)
velocityVector[0] = velocityX
velocityVector[1] = velocityY


local_shear_rate_X = FaceVariable(mesh=mesh, value = 0) # this term only acts 
on normal direction
local_shear_rate_Y = FaceVariable(mesh=mesh, value = 1) # linear velocity 
profile du/dz
local_shear_rate = FaceVariable(mesh=mesh, rank=1)
local_shear_rate[0] = local_shear_rate_X
local_shear_rate[1] = local_shear_rate_Y

Pe = 100. # some constants

#Setting up the diffusion part

D_Matrix = [[[0, 0],

 [0, 1/Pe]]]

#This is the new part to account for outflow fixed flux boundary (this does not 
set up correctly)
outflowflux = FaceVariable(mesh = mesh, value = phi.faceValue*velocityVector)
velocityVector.constrain(0.,mesh.facesRight)

# equation set up
eq = (TransientTerm() + UpwindConvectionTerm(coeff = velocityVector)\
+ UpwindConvectionTerm(coeff =local_shear_rate*(1-phi).faceValue)\
- (mesh.facesRight * outflowflux).divergence \
== DiffusionTerm (coeff=D_Matrix))


# set up the boundary conditions

phi.constrain(0.5, where = mesh.facesLeft) # fixed the inlet as 0.5 
concentration

phi.faceGrad.constrain(local_shear_rate*Pe*(1-phi.faceValue)*phi.faceValue, 
where = mesh.facesTop | mesh.facesBottom) #boundary conditions at the top and 
bottom


#phi.constrain(0, where = mesh.facesRight) #not working outlet boundary 
condition
#phi.faceGrad.constrain(velocityVector*phi.faceValue,where = mesh.facesRight) # 
problematic outlet boundary conditions

# Viewer set up

viewer = Matplotlib2

Re: Solving advection, segregation and diffusion equation

[Daniel Wheeler]

On Wed, Sep 14, 2016 at 7:52 PM, Guyer, Jonathan E. Dr. (Fed)
 wrote:
>
> I don't know offhand. With a Péclet number of 100, your problem is almost 
> completely hyperbolic, which FiPy (and cell-centered Finite Volume) isn't 
> very good at. Daniel knows more about this and may have some suggestions.
>
> You might consider adding a TransientTerm for artificial relaxation and 
> trying the VanLeerConvectionTerm.

Definitely use a TransientTerm and see the time step you can get away
with. Don't sweep, just use time steps. Using the VanLeerConvetionTerm
isn't necessary as the accuracy is not important in a steady state
problem.

-- 
Daniel Wheeler

___
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]

Re: Solving advection, segregation and diffusion equation

[Guyer, Jonathan E. Dr. (Fed)]

> On Sep 14, 2016, at 12:49 PM, Zhekai Deng  
> wrote:
> 
> So, in my previous one, I only have constraints on left (fixed 0.5 inlet 
> concentration), top and bottom ( the boundary conditions that I applied), but 
> not boundary condition on the right. I have implemented an additional one on 
> the right to account for the outlet:
> 
> phi.faceGrad.constrain(velocityVector*phi.faceValue,where = mesh.facesRight)
> 
> but it seems does not help very much.

I was surprised to see that several "corrections" I applied to your script 
didn't seem to change things very much. This problem is strangely insensitive 
to things that I think should matter.

In contrast, moving the definition of eq out of the sweep loop made a huge 
difference and I don't understand why.


> I am not sure what does "velocity is backwards someplace" mean. From the way 
> I defined it, the $u \hat{i} = \tilde{z}$ and $u 
> \hat{j} = 0$. Both are independent of solution variable. Would it be possible 
> to clarify what does "velocity is backwards" mean?

This was just idle speculation that maybe a velocity or shear points left when 
it should point right (up vs. down), etc. and that maybe that's what causes the 
rippling.


> After the discussion with my colleague, we think it is better to understand 
> this as a vector term that only act on the normal direction ($\hat{j}$) of 
> the flow domain. The reason for this is that the shear rate term is gradient 
> of velocity field which is a tensor. However, we are not entire sure the 
> functional dependence of segregation velocity on shear rate tensor. That's 
> why we usually call it as local shear rate (du/dz \hat{j} ), and I should 
> have been clear on this issue in the beginning.

Thanks for clarifying.


> In summary, it does seems that the initial sweeps give the right "trend" 
> toward the solution. However, I am still not sure about why it is not 
> convergent. Any suggestion on how to proceed to the next step ? I have 
> attached the newest version of the code.

I don't know offhand. With a Péclet number of 100, your problem is almost 
completely hyperbolic, which FiPy (and cell-centered Finite Volume) isn't very 
good at. Daniel knows more about this and may have some suggestions.

You might consider adding a TransientTerm for artificial relaxation and trying 
the VanLeerConvectionTerm.

- Jon

___
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]

Re: Solving advection, segregation and diffusion equation

[Zhekai Deng]

Thank you for the reply.  I try to answer any confusion in below.

1. I think what you have is fine. Explicitly, I'd write:

  phi.faceGrad.constrain(..., where=mesh.facesTop | mesh.facesBottom)

See http://www.ctcms.nist.gov/fipy/documentation/USAGE.html#
applying-spatially-varying-boundary-conditions

Thanks. It has been implemented in this way now.


2. For inlet/outlet conditions, please see the comment at

  http://www.ctcms.nist.gov/fipy/documentation/USAGE.html#
applying-outlet-or-inlet-boundary-conditions

If you impose constraints, then FiPy is inlet/outlet, not no-flux.

So, in my previous one, I only have constraints on left (fixed 0.5 inlet
concentration), top and bottom ( the boundary conditions that I applied),
but not boundary condition on the right. I have implemented an additional
one on the right to account for the outlet:

phi.faceGrad.constrain(velocityVector*phi.faceValue,where = mesh.facesRight)

but it seems does not help very much.

I've made some changes to your script that result in the initial sweeps
looking much more like your matlab solution, but it is still not
convergent. I wonder at this point if a velocity is backwards someplace.

Thank you for the change. Well, even after I added the right side outlet
boundary condition, the sweep solution still looks similar to what you
have. In the first initial sweeps, it looks like the matlab solution, but
it is not convergent. However, the residual seems does not explore as
quickly as before by adding the outlet condition.

I am not sure what does "velocity is backwards someplace" mean. From the
way I defined it, the $u \hat{i} = \tilde{z}$ and $u \hat{j} = 0$. Both are
independent of solution variable. Would it be possible to clarify what does
"velocity is backwards" mean?

I don't understand why $\tilde{\dot{\gamma}} = 1 \hat{\j}$.
Shouldn't it be $\tilde{\dot{\gamma}} = 1 \hat{\i}$, given that $\tilde{u}
= \tilde{z} \hat{\i}$?

Thanks for bring it up. After the discussion with my colleague, we think it
is better to understand this as a vector term that only act on the normal
direction ($\hat{j}$) of the flow domain. The reason for this is that the
shear rate term is gradient of velocity field which is a tensor. However,
we are not entire sure the functional dependence of segregation velocity on
shear rate tensor. That's why we usually call it as local shear rate (du/dz
\hat{j} ), and I should have been clear on this issue in the beginning.

Is, perhaps, the definition
  $\tilde{\dot{\gamma}} \equiv \frac{\partial \tilde{u}_x}{\partial
\tilde{z}} \hat{\j}$
where
  $\tilde{u} \equiv \tilde{u}_x \hat{\i} + \tilde{u}_y \hat{\j}$

Yes, this is the definition we are currently following.

In summary, it does seems that the initial sweeps give the right "trend"
toward the solution. However, I am still not sure about why it is not
convergent. Any suggestion on how to proceed to the next step ? I have
attached the newest version of the code.

Best,

Zhekai


On Tue, Sep 13, 2016 at 6:35 PM, Guyer, Jonathan E. Dr. (Fed) <
jonathan.gu...@nist.gov> wrote:

> Zhekai -
>
> 1. I think what you have is fine. Explicitly, I'd write:
>
>   phi.faceGrad.constrain(..., where=mesh.facesTop | mesh.facesBottom)
>
> See http://www.ctcms.nist.gov/fipy/documentation/USAGE.html#
> applying-spatially-varying-boundary-conditions
>
>
> 2. For inlet/outlet conditions, please see the comment at
>
>   http://www.ctcms.nist.gov/fipy/documentation/USAGE.html#
> applying-outlet-or-inlet-boundary-conditions
>
> If you impose constraints, then FiPy is inlet/outlet, not no-flux.
>
>
> I've made some changes to your script that result in the initial sweeps
> looking much more like your matlab solution, but it is still not
> convergent. I wonder at this point if a velocity is backwards someplace.
>
> I don't understand why $\tilde{\dot{\gamma}} = 1 \hat{\j}$.
> Shouldn't it be $\tilde{\dot{\gamma}} = 1 \hat{\i}$, given that $\tilde{u}
> = \tilde{z} \hat{\i}$?
>
> Is, perhaps, the definition
>   $\tilde{\dot{\gamma}} \equiv \frac{\partial \tilde{u}_x}{\partial
> \tilde{z}} \hat{\j}$
> where
>   $\tilde{u} \equiv \tilde{u}_x \hat{\i} + \tilde{u}_y \hat{\j}$
>
> (I *really* strongly recommend not writing vectors like scalars)
>
> - Jon
>
> > On Sep 13, 2016, at 5:45 PM, Zhekai Deng  northwestern.edu> wrote:
> >
> > Hello All,
> >
> > I am trying to use Fipy to solve advection, segregation and diffusion
> equation in 2D. I have attached the pdf to explain in details about my
> equation, and corresponding fipy code that I have problem with it.
> >
> > The solution I got has "kind of similar shape" to the same equation that
> I solved using MATLAB(the code is also attached, I have verified that
> MATLAB is solving correctly). However, the code I have in Fipy does not
> converge, and the residual jump between large and small number.
> >
> > Also,I have a few things that I am not entire sure.
> >
> > 1. In 2D, how could we specify the which direction of the
> "faceGrad.constrain" t

Re: Solving advection, segregation and diffusion equation

[Guyer, Jonathan E. Dr. (Fed)]

Zhekai -

1. I think what you have is fine. Explicitly, I'd write:

  phi.faceGrad.constrain(..., where=mesh.facesTop | mesh.facesBottom)

See 
http://www.ctcms.nist.gov/fipy/documentation/USAGE.html#applying-spatially-varying-boundary-conditions


2. For inlet/outlet conditions, please see the comment at 

  
http://www.ctcms.nist.gov/fipy/documentation/USAGE.html#applying-outlet-or-inlet-boundary-conditions

If you impose constraints, then FiPy is inlet/outlet, not no-flux.


I've made some changes to your script that result in the initial sweeps looking 
much more like your matlab solution, but it is still not convergent. I wonder 
at this point if a velocity is backwards someplace.

I don't understand why $\tilde{\dot{\gamma}} = 1 \hat{\j}$. 
Shouldn't it be $\tilde{\dot{\gamma}} = 1 \hat{\i}$, given that $\tilde{u} = 
\tilde{z} \hat{\i}$?

Is, perhaps, the definition 
  $\tilde{\dot{\gamma}} \equiv \frac{\partial \tilde{u}_x}{\partial \tilde{z}} 
\hat{\j}$
where 
  $\tilde{u} \equiv \tilde{u}_x \hat{\i} + \tilde{u}_y \hat{\j}$

(I *really* strongly recommend not writing vectors like scalars)

- Jon

> On Sep 13, 2016, at 5:45 PM, Zhekai Deng  
> wrote:
> 
> Hello All,
> 
> I am trying to use Fipy to solve advection, segregation and diffusion 
> equation in 2D. I have attached the pdf to explain in details about my 
> equation, and corresponding fipy code that I have problem with it.
> 
> The solution I got has "kind of similar shape" to the same equation that I 
> solved using MATLAB(the code is also attached, I have verified that MATLAB is 
> solving correctly). However, the code I have in Fipy does not converge, and 
> the residual jump between large and small number.
> 
> Also,I have a few things that I am not entire sure.
> 
> 1. In 2D, how could we specify the which direction of the 
> "faceGrad.constrain" that we want to constrain ?
> 2. In my problem specifically,  I would like to implement convection inlet 
> flux and convection outlet flux. This seems is against the Fipy default no 
> flux condition.  Any suggestion to on how to implement it ? 
> 
> Best,
> 
> Zhekai
> ___
> fipy mailing list
> fipy@nist.gov
> http://www.ctcms.nist.gov/fipy
>  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
from fipy import *

# mesh set up

nx = 300.

nz = 300.

dx = 1/nx

dz = 1/nz



mesh = Grid2D(dx=dx, dy=dz, nx=nx, ny=nz)

phi = CellVariable(name = "solution variable",mesh = mesh, value = 0.5)



# the follwing set up the velocity and shear rate in vector form.

# for those direction that I don't need it, I set its value as zero.

velocityX = mesh.faceCenters[1] # linear velocity profile

velocityY = FaceVariable(mesh=mesh, value = 0) # zero velocity in y direction

velocityVector = FaceVariable(mesh=mesh, rank=1)

velocityVector[0] = velocityX

velocityVector[1] = velocityY



shear_rate_X = FaceVariable(mesh=mesh, value = 0) # linear velocity profile

shear_rate_Y = FaceVariable(mesh=mesh, value = 1) # zero velocity in y direction

shear_rate = FaceVariable(mesh=mesh, rank=1)

shear_rate[0] = shear_rate_X

shear_rate[1] = shear_rate_Y





Pe = 100. # some constants



#Setting up the diffusion part

D_Matrix = [[[0, 0],

 [0, 1/Pe]]]



eq = (ExponentialConvectionTerm(coeff = velocityVector) + ExponentialConvectionTerm(coeff = shear_rate*(1-phi).faceValue) == DiffusionTerm (coeff=D_Matrix))



# mask represents where the boundary condition is applied

mask = mesh.facesTop | mesh.facesBottom 

phi.constrain(0.5, where = mesh.facesLeft) # fixed the inlet as 0.5 concentration

phi.faceGrad.constrain(shear_rate*Pe*(1-phi.faceValue)*phi.faceValue,mask) #boundary conditions at the top and bottom (Not sure whether it is correct or not)



# Viewer set up

viewer = Matplotlib2DViewer(vars=phi,datamin=0.,datamax = 1., title="final solution")



# Begin the calculation of the pde

res = 1e+10

while res >1e-6: #wait until it converges

res = eq.sweep(var = phi)

viewer.plot()

print res

#I implmented both advection and segregation term using ExponentialConvectionTerm. 

#However, I would want there is a inlet flux on the left, and outlet flux on the right due to only advection.

#Any suggestions on how to implement is very appreciated! 

___
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]

Solving advection, segregation and diffusion equation

[Zhekai Deng]

Hello All,

I am trying to use Fipy to solve advection, segregation and diffusion
equation in 2D. I have attached the pdf to explain in details about my
equation, and corresponding fipy code that I have problem with it.

The solution I got has "kind of similar shape" to the same equation that I
solved using MATLAB(the code is also attached, I have verified that MATLAB
is solving correctly). However, the code I have in Fipy does not converge,
and the residual jump between large and small number.

Also,I have a few things that I am not entire sure.

1. In 2D, how could we specify the which direction of the
"faceGrad.constrain" that we want to constrain ?
2. In my problem specifically,  I would like to implement convection inlet
flux and convection outlet flux. This seems is against the Fipy default no
flux condition.  Any suggestion to on how to implement it ?

Best,

Zhekai
from fipy import *
# mesh set up
nx = 300.
nz = 300.
dx = 1/nx
dz = 1/nz

mesh = Grid2D(dx=dx, dy=dz, nx=nx, ny=nz)
phi = CellVariable(name = "solution variable",mesh = mesh, value = 0.5)

# the follwing set up the velocity and shear rate in vector form.
# for those direction that I don't need it, I set its value as zero.
velocityX = mesh.faceCenters[1] # linear velocity profile
velocityY = FaceVariable(mesh=mesh, value = 0) # zero velocity in y direction
velocityVector = FaceVariable(mesh=mesh, rank=1)
velocityVector[0] = velocityX
velocityVector[1] = velocityY

shear_rate_X = FaceVariable(mesh=mesh, value = 0) # linear velocity profile
shear_rate_Y = FaceVariable(mesh=mesh, value = 1) # zero velocity in y direction
shear_rate = FaceVariable(mesh=mesh, rank=1)
shear_rate[0] = shear_rate_X
shear_rate[1] = shear_rate_Y


Pe = 100. # some constants

#Setting up the diffusion part
Diffusion_X = FaceVariable(mesh=mesh, value = 0)
Diffusion_Y = FaceVariable(mesh=mesh, value = 1/Pe)
D_Vector = FaceVariable(mesh=mesh, rank=1)
D_Vector[0] = Diffusion_X
D_Vector[1] = Diffusion_Y

# mask represents where the boundary condition is applied
mask = mesh.facesTop | mesh.facesBottom 
phi.constrain(0.5, where = mesh.facesLeft) # fixed the inlet as 0.5 
concentration
phi.faceGrad.constrain([shear_rate*Pe*(1-phi.faceValue)*phi.faceValue],mask) 
#boundary conditions at the top and bottom (Not sure whether it is correct or 
not)

# Viewer set up
viewer = Matplotlib2DViewer(vars=phi,datamin=0.,datamax = 1., title="final 
solution")

# Begin the calculation of the pde
res = 1e+10
while res >1e-6: #wait until it converges
eq = (ExponentialConvectionTerm(coeff = velocityVector) + 
ExponentialConvectionTerm(coeff = shear_rate*(1-phi).faceValue) == 
DiffusionTerm (coeff=D_Vector))
res = eq.sweep(var = phi)
viewer.plot()
print res
#I implmented both advection and segregation term using 
ExponentialConvectionTerm. 
#However, I would want there is a inlet flux on the left, and outlet flux on 
the right due to only advection.
#Any suggestions on how to implement is very appreciated! 


matlab_demo.m
Description: application/vnd.wolfram.mathematica.package


equation_explain.pdf
Description: Adobe PDF document
___
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]