Re: ternary electrolyte
;>>> >>>>>>> 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 >>>>>>>> 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 >>>>>>>>> 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 >>>>>>>>>> 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 equa
Re: ternary electrolyte
nks 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 >>>>>>>>> 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 >>>>>>>>>> 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=http-3A__bit.ly_2wJDXob-3D-3D=DwIGaQ=nE__W8dFE-shTxStwXtp0A=GFdsMpHIw-aduHn9jOIpyQ=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c=7BJyGej3bHZQwW3KYHek8AbNpqq32E63p3rH-Vvq9UM= >>>>>>>>>>
Re: ternary electrolyte
>>>>> 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 >>>>>> 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 >>>>>>> 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 >>>>>>>> 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 >>>>>>>&g
Re: ternary electrolyte
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 >>>>> mailto:fipy@nist.gov>> 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 >>>>> <mailto:s...@msu.edu>> 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 >>>>>>> mailto:fipy@nist.gov>> 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 n
Re: ternary electrolyte
;> 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 >>>>>> 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 >>>>>>> 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 http://bit.ly/2wJDXob== >>>>>>> >>>>>>> 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: >>>>>>> http://bit.ly/2wOzxfK== >>>>>>> >>>>>>> Any insights or suggestions would be welcome. >>>>>>> >>>>>>> Cheers, >>>>>>> Scott☘ >>>>>>> >>>>>>> >>>>>>> ___ >>>>>>> fipy mailing list >>>>>>> fipy@nist.gov >>>>>>> http://www.ctcms.nist.gov/fipy== >>>>>>> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy== >>>>>>> ] >>>>>> >>>>>> >>>>>> ___ >>>>>> fipy mailing list >>>>>> fipy@nist.gov >>>>>> http://www.ctcms.nist.gov/fipy== >>>>>> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy==] >>>>>> >>>>> >>>>> >>>>> ___ >>>>> fipy mailing list >>>>> fipy@nist.gov >>>>> http://www.ctcms.nist.gov/fipy== >>>>> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy== ] >>>> >>>> >>>> ___ >>>> fipy mailing list >>>> fipy@nist.gov >>>> http://www.ctcms.nist.gov/fipy== >>>> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy== ] >>>> >>> >> >> >> ___ >> fipy mailing list >> fipy@nist.gov >> http://www.ctcms.nist.gov/fipy= >> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy= ] >> > ___ fipy mailing list fipy@nist.gov http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
Re: ternary electrolyte
; 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 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=http-3A__bit.ly_2wJDXob-3D=DwIGaQ=nE__W8dFE-shTxStwXtp0A=GFdsMpHIw-aduHn9jOIpyQ=KNN_W5o3-EoZzqSQCJS2Ra1FTVyZod7m5mC_mKj6USY=DgoS_BGUnGcKKylyNkhKcWC11wLatFOb5VRN9YYOuwU= >>>>>> >>>>>> >>>>>> 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 diffu
Re: ternary electrolyte
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 > 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 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-257C1687dc7e080545a63f0708d6eeae58e7-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C1-257C636958829930551626-26amp-3Bsdata-3DcSBGKKJG0-252Ft-252FMtTZbnTwenKAVsKPHa-252F1C8NgupkRI0E-253D-26amp-3Breserved-3D0=DwIGaQ=nE__W8dFE-shTxStwXtp0A=GFdsMpHIw-aduHn9jOIpyQ=MiNELfmZGpKXQZs2VXW1BGZb0xHBiSGDtwc09TBIN_4=mskW-yZCCtkB6VJ1e2iNv6InRpjYLWWU0Xklumrtr0o= >> >> >> 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-257C1687dc7e080545a63f0708d6eeae58e7-257C2ab5d82fd8fa4797a93e054655c61dec-257C1-257C1-257C636958829930
Re: ternary electrolyte
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 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://gcc01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fbit.ly%2F2wJDXobdata=02%7C01%7Cjonathan.guyer%40nist.gov%7C1687dc7e080545a63f0708d6eeae58e7%7C2ab5d82fd8fa4797a93e054655c61dec%7C1%7C1%7C636958829930551626sdata=cSBGKKJG0%2Ft%2FMtTZbnTwenKAVsKPHa%2F1C8NgupkRI0E%3Dreserved=0 > > 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://gcc01.safelinks.protection.outlook.com/?url=http%3A%2F%2Fbit.ly%2F2wOzxfKdata=02%7C01%7Cjonathan.guyer%40nist.gov%7C1687dc7e080545a63f0708d6eeae58e7%7C2ab5d82fd8fa4797a93e054655c61dec%7C1%7C1%7C636958829930551626sdata=CkRvDWJQl%2FjHRTwtxd6y1Krs3VuALc2bAOXKJYsdmEk%3Dreserved=0 > > Any insights or suggestions would be welcome. > > Cheers, > Scott☘ > > > ___ > fipy mailing list > fipy@nist.gov > http://www.ctcms.nist.gov/fipy > [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ] ___ fipy mailing list fipy@nist.gov http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
ternary electrolyte
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 http://bit.ly/2wJDXob 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: http://bit.ly/2wOzxfK Any insights or suggestions would be welcome. Cheers, Scott☘ ___ fipy mailing list fipy@nist.gov http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
