1. As before, there is no flux condition on species 3 2. Unless your transport is at the speed of light, nonzero flux will not alter Poisson's equation. That's probably a red herring, though, for now.
More germane, the notebook at http://bit.ly/2wOzxfK does not solve for me, with a singular matrix error. Replacing C2.faceGrad.divergence with DiffusionTerm(coeff=-z1/z2, var=C1) gives a converged solution, but again, Phi = 0. Phi = 0, C1 = C2, with grad(C1) = grad(C2) = const is a solution to this set of equations. I still see nothing in the conditions you have provided that would induce Phi to curve. I have tried both SciPy and Pysparse and I don't get the results you show for the binary. Can you specify what environment you're running on? > On Jun 14, 2019, at 1:32 PM, Scott Calabrese Barton <s...@msu.edu> wrote: > > Hi, two answers: > > 1. Phi=0 is not a solution because it does not satisfy the flux condition on > species 3. > > 2. If nonuniform charge or nonuniform permittivity are my only two choices, > I’d go with nonuniform permittivity, due to varying concentrations. However, > this is non-equilibrium system so I would suggest that nonzero flux can also > give rise to potential gradients. This may be more clear in the binary system > http://bit.ly/2wOzxfK. > > A physical manifestation of this system is copper deposition across a > boundary layer, say from copper sulfate with a supporting electrolyte such as > sulfuric acid. In this case the three species might be copper, protons, and > sulfate ions. > > ☘ > >> On Jun 14, 2019, at 11:26 AM, Guyer, Jonathan E. Dr. (Fed) via fipy >> <fipy@nist.gov> wrote: >> >> Phi = 0 is a solution. Are you suggesting that this set of equations has >> more than one solution? >> >> I realize you are not solving Poisson's equation. My question was to >> consider the broader physics of the problem. Physically, non-uniform >> electric field means either that charge is non-uniform, which is not in your >> problem, or permittivity is non-uniform, which is not in your problem. So, >> what's in your problem that resolves that physical inconsistency?. >> >>> On Jun 14, 2019, at 10:53 AM, Scott Calabrese Barton <s...@msu.edu> wrote: >>> >>> To answer your questions first, the non-zero gradient in Phi counteracts >>> the gradients in C2 and C3 such that the migration and diffusion flux >>> cancel and overall flux is zero for those species. Regarding the curvature, >>> you’re thinking of the Poisson equation that is only valid in this context >>> when concentrations are constant. The variations in concentration give >>> rise to nonlinear variations in Phi. >>> >>> Perhaps I should point out that we can get the numerical solution using >>> other solvers such as solve_bvp from scipy, and an analytical solution is >>> also available (that’s the source of the plot). So the problem is well >>> defined. I’ll be happy to share those calculations if needed. I’m using >>> this problem as an example to learn more about fipy. >>> >>>> On Jun 14, 2019, at 9:10 AM, Guyer, Jonathan E. Dr. (Fed) via fipy >>>> <fipy@nist.gov> wrote: >>>> >>>> I had the same initial thought, but you're not adding new information by >>>> doing that. You can see this by looking at the matrix stencil that results >>>> from these equations: >>>> >>>> X X C1 >>>> X C2 >>>> XX C3 = 0 >>>> XXX Phi >>>> >>>> That bottom right sector of the matrix is empty and so the matrix solvers >>>> cannot find a solution because you're now trying to solve for four >>>> unknowns with what are really only three equations. >>>> >>>> Things to consider: >>>> >>>> In the solution you posted from West, Phi does not have zero gradient at >>>> either end of the domain. Why? What boundary condition is driving that? >>>> >>>> Phi is curved, which implies that there is net charge in the system, but >>>> you assert zero charge throughout. How does West reconcile that >>>> inconsistency? >>>> >>>>> On Jun 13, 2019, at 11:25 AM, Scott Calabrese Barton <s...@msu.edu> wrote: >>>>> >>>>> That’s true, so we can make C3 explicit (full code at >>>>> https://urldefense.proofpoint.com/v2/url?u=http-3A__bit.ly_2F7UQgR-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=3mcnknSVK1Ij82U2PK6itcevIOeNgXVMqLhUYYihMqc&e= >>>>> ): >>>>> >>>>> # 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 >>>>>> <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 <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 >>>>>>>> <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 need some form of Poisson's equation. >>>>>>>> >>>>>>>> - Jon >>>>>>>> >>>>>>>>> On Jun 11, 2019, at 4:48 PM, Scott Calabrese Barton <s...@msu.edu> >>>>>>>>> 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&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=7BJyGej3bHZQwW3KYHek8AbNpqq32E63p3rH-Vvq9UM&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=http-3A__bit.ly_2wOzxfK-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=m59IS6E5x3ikFID0-LV2HY3XwhWMEDBSf7p0NOFfTSI&e= >>>>>>>>> >>>>>>>>> >>>>>>>>> Any insights or suggestions would be welcome. >>>>>>>>> >>>>>>>>> Cheers, >>>>>>>>> Scott☘ >>>>>>>>> >>>>>>>>> >>>>>>>>> _______________________________________________ >>>>>>>>> fipy mailing list >>>>>>>>> fipy@nist.gov >>>>>>>>> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=ReSWGMkuYcc6b3YlxVB_e86OwfNxqVxp12I2NbAXyyM&e= >>>>>>>>> >>>>>>>>> [ NIST internal ONLY: >>>>>>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__email.nist.gov_mailman_listinfo_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=yCpx4zcj6bXxIF1fxMvVdOKlrt0cB2KI4MtVr_GxD8E&e= >>>>>>>>> ] >>>>>>>> >>>>>>>> >>>>>>>> _______________________________________________ >>>>>>>> fipy mailing list >>>>>>>> fipy@nist.gov >>>>>>>> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=ReSWGMkuYcc6b3YlxVB_e86OwfNxqVxp12I2NbAXyyM&e= >>>>>>>> >>>>>>>> [ NIST internal ONLY: >>>>>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__email.nist.gov_mailman_listinfo_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=yCpx4zcj6bXxIF1fxMvVdOKlrt0cB2KI4MtVr_GxD8E&e= >>>>>>>> ] >>>>>>>> >>>>>>> >>>>>>> >>>>>>> _______________________________________________ >>>>>>> fipy mailing list >>>>>>> fipy@nist.gov >>>>>>> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=ReSWGMkuYcc6b3YlxVB_e86OwfNxqVxp12I2NbAXyyM&e= >>>>>>> >>>>>>> [ NIST internal ONLY: >>>>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__email.nist.gov_mailman_listinfo_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=yCpx4zcj6bXxIF1fxMvVdOKlrt0cB2KI4MtVr_GxD8E&e= >>>>>>> ] >>>>>> >>>>>> >>>>>> _______________________________________________ >>>>>> fipy mailing list >>>>>> fipy@nist.gov >>>>>> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=ReSWGMkuYcc6b3YlxVB_e86OwfNxqVxp12I2NbAXyyM&e= >>>>>> >>>>>> [ NIST internal ONLY: >>>>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__email.nist.gov_mailman_listinfo_fipy-3D-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=yCpx4zcj6bXxIF1fxMvVdOKlrt0cB2KI4MtVr_GxD8E&e= >>>>>> ] >>>>>> >>>>> >>>> >>>> >>>> _______________________________________________ >>>> fipy mailing list >>>> fipy@nist.gov >>>> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=33cYoWvm54o3WvNWYWRvqnGY0sRat3jfaXk8q_X5w9I&e= >>>> >>>> [ NIST internal ONLY: >>>> https://urldefense.proofpoint.com/v2/url?u=https-3A__email.nist.gov_mailman_listinfo_fipy-3D&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=XRcoJVJWV1WhFkofW6wT0XWCpKQQCqBAbqHQmFaT-Zg&e= >>>> ] >>>> >>> >> >> >> _______________________________________________ >> fipy mailing list >> fipy@nist.gov >> https://urldefense.proofpoint.com/v2/url?u=http-3A__www.ctcms.nist.gov_fipy&d=DwIGaQ&c=nE__W8dFE-shTxStwXtp0A&r=GFdsMpHIw-aduHn9jOIpyQ&m=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=eqC4_t0f-3FwcRlXJXk2ZHnlT5ywI795EoHNf7OTS1o&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=NXf47PawWEovDQCuiqrE25DZuppCDLEVN9L9lPfcl3c&s=V8ZHNRbUJRgcU2d8gY0CNkaGSM-rSRB_-VfVsP_CPF8&e= >> ] _______________________________________________ fipy mailing list fipy@nist.gov http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
