That’s true, so we can make C3 explicit (full code at http://bit.ly/2F7UQgR):
# Variables
C1= CellVariable(mesh = mesh, value = 1.0)
C2= CellVariable(mesh = mesh, value = 4.0/3.0)
Phi = CellVariable(mesh = mesh, value = 0)
C30 = -1/z3 * (z1*1.0 + z2*C20) # Electroneutrality
C3= CellVariable(mesh = mesh, value = C30)
# Boundary Conditions
C1.constrain(1.0, mesh.facesRight)
C2.constrain(C20, mesh.facesRight)
C3.constrain(C30, mesh.facesRight)
Phi.constrain(-0.2, mesh.facesRight)
C1.constrain(0.0, mesh.facesLeft)
# Governing Equations
Eq1 = + DiffusionTerm(coeff=z1*C1, var=Phi) + DiffusionTerm(coeff=1.0, var=C1)
Eq2 = + DiffusionTerm(coeff=z2*C2, var=Phi) + DiffusionTerm(coeff=1.0, var=C2)
Eq3 = + DiffusionTerm(coeff=z3*C3, var=Phi) + DiffusionTerm(coeff=1.0, var=C3)
Eq4 = ImplicitSourceTerm(coeff=z1, var=C1) + \
ImplicitSourceTerm(coeff=z2, var=C2) + \
ImplicitSourceTerm(coeff=z3, var=C3)
Eqns = Eq1 & Eq2 & Eq3 & Eq4
However, this set of equations does not converge. I’m fairly certain that I am
mismanaging the boundary conditions.
> On Jun 12, 2019, at 5:02 PM, Guyer, Jonathan E. Dr. (Fed) via fipy
> <[email protected]> wrote:
>
> There's no zero-flux BC on C3 because there's no such thing as C3 as far as
> FiPy is concerned. The boundary conditions on C3 are entirely determined by
> the boundary conditions on C1 and C2 (Dirichlet on the right and Robin on the
> left). The solution obtained satisfies the equations and boundary conditions
> given to FiPy.
>
>> On Jun 12, 2019, at 2:19 PM, Scott Calabrese Barton <[email protected]> wrote:
>>
>> Jon, thanks for catching that typo. I wish it were the problem!
>>
>> The solution with grad phi=0 does not satisfy the zero-flux boundary
>> condition on species 3. In the absence of a potential gradient, that
>> species would need to also have zero gradient (like species 2 in the current
>> solution).
>>
>> So the question becomes how to enforce that zero-flux BC. Any suggestions?
>> I could drop the conservation equation and replace it with something that
>> enforces zero-flux everywhere, but I’d prefer to use a boundary condition as
>> that would be more applicable in higher dimensions.
>>
>> ☘
>>
>>
>>> On Jun 12, 2019, at 1:36 PM, Guyer, Jonathan E. Dr. (Fed) via fipy
>>> <[email protected]> wrote:
>>>
>>> Scott -
>>>
>>> Eq2 appears to have a typo:
>>>
>>> -Eq2 = DiffusionTerm(coeff=z2*C1, var=Phi) + DiffusionTerm(coeff=1.0,
>>> var=C2)
>>> +Eq2 = DiffusionTerm(coeff=z2*C2, var=Phi) + DiffusionTerm(coeff=1.0,
>>> var=C2)
>>>
>>> Eq3 has (much) better coupling if written like this:
>>>
>>> -Eq3 = DiffusionTerm(coeff=z3*C3, var=Phi) + (C3.faceGrad).divergence
>>> +Eq3 = DiffusionTerm(coeff=z3*C3, var=Phi) - DiffusionTerm(coeff=-z1/z3,
>>> var=C1) - DiffusionTerm(coeff=-z2/z3, var=C2)
>>>
>>> [Note that this way you don't need to evaluate .faceGrad at the boundaries,
>>> so your issue #650 doesn't come up. We don't upwind at the boundaries
>>> because it's complicated and expensive to do right and it's unusual to
>>> care.]
>>>
>>> With these changes, the solution converges in one sweep, but does not agree
>>> with the Analytical Solution you give from West.
>>>
>>> I don't understand a couple of things though. As posed, what governs \Phi?
>>> I get \Phi = 0 and see no reason why it would be otherwise. Typically, I
>>> would expect Poisson's equation to be obeyed. Since you assert charge
>>> electroneutrality, this reduces to Laplace's equation, but in that case,
>>> the slope of \Phi should be constant (unless the permittivity is
>>> non-uniform, but you don't have a permittivity). In the Analytical plot,
>>> \Phi has definite curvature.
>>>
>>> My guess is you don't need Eq3 (it'll get taken care of automatically) and
>>> you need some form of Poisson's equation.
>>>
>>> - Jon
>>>
>>>> On Jun 11, 2019, at 4:48 PM, Scott Calabrese Barton <[email protected]> wrote:
>>>>
>>>> I’ve had some trouble using fipy for the 1D, steady-state ternary
>>>> electrolyte problem, which involves coupled, nonlinear migration and
>>>> diffusion.
>>>>
>>>> A text version of the model is listed below. A more complete development
>>>> is at
>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__gcc01.safelinks.protection.outlook.com_-3Furl-3Dhttp-253A-252F-252Fbit.ly-252F2wJDXob-26amp-3Bdata-3D02-257C01-257Cjonathan.guyer-2540nist.gov-257C76b3f81e57154ae5f33f08d6ef62b1ed-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C1-257C636959604508701127-26amp-3Bsdata-3DIVWC4t-252F3V908aQtRyOmUCUHFVxMn96-252BZfyE457gURv0-253D-26amp-3Breserved-3D0-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=8qFl8yyF0LbTlYEFAOX0uymiY7LalYSuAcmcV0ozSNg&e=
>>>>
>>>>
>>>> from fipy import *
>>>>
>>>> # Parameters
>>>> nx=50
>>>> C20=4.0/3.0
>>>> z1, z2, z3= 2.0, -2.0, 1.0
>>>>
>>>> # Mesh
>>>> mesh = Grid1D(nx=nx, dx=1.0/nx)
>>>>
>>>> # Variables
>>>> C1= CellVariable(mesh = mesh, value = 1.0)
>>>> C2= CellVariable(mesh = mesh, value = 4.0/3.0)
>>>> Phi = CellVariable(mesh = mesh, value = -0.5)
>>>>
>>>> # Electroneutrality
>>>> C3 = -1/z3 * (z1*C1 + z2*C2)
>>>>
>>>> # Boundary Conditions
>>>> C1.constrain(1.0, mesh.facesRight)
>>>> C1.constrain(0.0, mesh.facesLeft)
>>>> C2.constrain(C20, mesh.facesRight)
>>>> Phi.constrain(0, mesh.facesRight)
>>>>
>>>> # Governing Equations
>>>> Eq1 = DiffusionTerm(coeff=z1*C1, var=Phi) + DiffusionTerm(coeff=1.0,
>>>> var=C1)
>>>> Eq2 = DiffusionTerm(coeff=z2*C1, var=Phi) + DiffusionTerm(coeff=1.0,
>>>> var=C2)
>>>> Eq3 = DiffusionTerm(coeff=z3*C3, var=Phi) + (C3.faceGrad).divergence
>>>>
>>>> Eqns = Eq1 & Eq2 & Eq3
>>>>
>>>> # Solution
>>>> res = 1e+10
>>>> restol= 1e-3
>>>>
>>>> while res > restol:
>>>> res = Eqns.sweep()
>>>> print(res)
>>>>
>>>>
>>>> However, this system does not converge. I suspect the issue is with
>>>> boundary conditions. I’ve tried zeroing out the diffusion coefficients at
>>>> the boundary, but that does not seem to make a difference here. I am able
>>>> to use this approach for the analogous binary problem:
>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__gcc01.safelinks.protection.outlook.com_-3Furl-3Dhttp-253A-252F-252Fbit.ly-252F2wOzxfK-26amp-3Bdata-3D02-257C01-257Cjonathan.guyer-2540nist.gov-257C76b3f81e57154ae5f33f08d6ef62b1ed-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C1-257C636959604508701127-26amp-3Bsdata-3DfPKk1duJaGEmdPa5zgSW4KPZmLzjp-252FRZ4mPi5xNMAcE-253D-26amp-3Breserved-3D0-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=9y0e0cL4FxuXNO8wKVX1mBA2uowCpCPyp3wfg5QcrbQ&e=
>>>>
>>>>
>>>> Any insights or suggestions would be welcome.
>>>>
>>>> Cheers,
>>>> Scott☘
>>>>
>>>>
>>>> _______________________________________________
>>>> fipy mailing list
>>>> [email protected]
>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__gcc01.safelinks.protection.outlook.com_-3Furl-3Dhttp-253A-252F-252Fwww.ctcms.nist.gov-252Ffipy-26amp-3Bdata-3D02-257C01-257Cjonathan.guyer-2540nist.gov-257C76b3f81e57154ae5f33f08d6ef62b1ed-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C0-257C636959604508701127-26amp-3Bsdata-3DJ7W3yaiabC4sNjcBUwB-252BEkTFZ-252Bd9UwBneyURQwikVGs-253D-26amp-3Breserved-3D0-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=r6OiuU-YfGLO-Eqgaw5BQWLseUqSMpLEl3UFSbkxCYw&e=
>>>>
>>>> [ NIST internal ONLY:
>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__gcc01.safelinks.protection.outlook.com_-3Furl-3Dhttps-253A-252F-252Femail.nist.gov-252Fmailman-252Flistinfo-252Ffipy-26amp-3Bdata-3D02-257C01-257Cjonathan.guyer-2540nist.gov-257C76b3f81e57154ae5f33f08d6ef62b1ed-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C0-257C636959604508701127-26amp-3Bsdata-3D7RzkhG-252FAIf-252BhmhFTfZuaDXqaPZWu9FXMXrslkgxvf3c-253D-26amp-3Breserved-3D0-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=fYyJlClOGv3gPUssX-7-NOiXZQJh-xUDvm8sExNXOGo&e=
>>>> ]
>>>
>>>
>>> _______________________________________________
>>> fipy mailing list
>>> [email protected]
>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__gcc01.safelinks.protection.outlook.com_-3Furl-3Dhttp-253A-252F-252Fwww.ctcms.nist.gov-252Ffipy-26amp-3Bdata-3D02-257C01-257Cjonathan.guyer-2540nist.gov-257C76b3f81e57154ae5f33f08d6ef62b1ed-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C0-257C636959604508701127-26amp-3Bsdata-3DJ7W3yaiabC4sNjcBUwB-252BEkTFZ-252Bd9UwBneyURQwikVGs-253D-26amp-3Breserved-3D0-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=r6OiuU-YfGLO-Eqgaw5BQWLseUqSMpLEl3UFSbkxCYw&e=
>>>
>>> [ NIST internal ONLY:
>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__gcc01.safelinks.protection.outlook.com_-3Furl-3Dhttps-253A-252F-252Femail.nist.gov-252Fmailman-252Flistinfo-252Ffipy-26amp-3Bdata-3D02-257C01-257Cjonathan.guyer-2540nist.gov-257C76b3f81e57154ae5f33f08d6ef62b1ed-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C0-257C636959604508701127-26amp-3Bsdata-3D7RzkhG-252FAIf-252BhmhFTfZuaDXqaPZWu9FXMXrslkgxvf3c-253D-26amp-3Breserved-3D0-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=fYyJlClOGv3gPUssX-7-NOiXZQJh-xUDvm8sExNXOGo&e=
>>> ]
>>>
>>
>>
>> _______________________________________________
>> fipy mailing list
>> [email protected]
>> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=QzZR5WAg-R5UuBOWw99RX4R7O3Ih6h0u1KkM4OGyOig&e=
>>
>> [ NIST internal ONLY:
>> https://urldefense.proofpoint.com/v2/url?u=https-3A__email.nist.gov_mailman_listinfo_fipy&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=klbQX_KDEfh5xGDVaXPT4Z_aevADif1RfbzdE1SPvoo&e=
>> ]
>
>
> _______________________________________________
> fipy mailing list
> [email protected]
> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=QzZR5WAg-R5UuBOWw99RX4R7O3Ih6h0u1KkM4OGyOig&e=
>
> [ NIST internal ONLY:
> https://urldefense.proofpoint.com/v2/url?u=https-3A__email.nist.gov_mailman_listinfo_fipy&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=KWl6oWKQ_cDGcu9bvRcy5gJc9Qd9TX9fm8VqHThYIRg&s=klbQX_KDEfh5xGDVaXPT4Z_aevADif1RfbzdE1SPvoo&e=
> ]
>
_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
