I'm not sure why you conclude from that example that var.faceGrad.constrain(0, where=m.exteriorFaces)
only applies to the advection component. Where does it say that? In the code you sent, with eq3 = TransientTerm() == DiffusionTerm(coeff=D)+ convection(coeff=faceVelocity) the flux is D \nabla \phi + \phi b \tanh(alpha*r) (pos-den)/r By constraining var.faceGrad to zero, you are saying that the boundary flux is \phi b \tanh(alpha*r) (pos-den)/r This is not a conservative boundary condition. > On Aug 24, 2016, at 8:59 PM, Yun Tao <yun...@ucdavis.edu> wrote: > > Hi Dan, > > Wow. I did not expect that at all, thanks for that important information! I > just tried searching for ways to implement the Robin BCs in FiPy. This site > example > <http://www.ctcms.nist.gov/fipy/examples/convection/generated/examples.convection.robin.html> > appears to show that the var.faceGrad.constrain(0, where=m.exteriorFaces) > command I mentioned earlier handles only the advection component, while, I > presume, the diffusion component is subject to Neumann BCs by default. If > that is so, then it's even more puzzling why having that line in my code > blows up the solution, yet the solution seems well-behaved under just the > Neumann BCs. Am I missing something else here? > > Thanks, > Yun > > On Wed, Aug 24, 2016 at 8:24 PM, Daniel Farrell <boyfarr...@gmail.com> wrote: > Hello Yun, > > I just briefly looked at the code. Seems like you are solving something like > an advection diffusion problem. This might not help but just I case ... > > If you want a closed boundary you need to apply Robin boundary conditions > because by definition the flux contains two components: one related to the > diffusion process and one the advection process. > > For example, http://scicomp.stackexchange.com/a/10576/3691 > > Dan > > On 24 Aug 2016, at 22:56, Yun Tao <yun...@ucdavis.edu> wrote: > >> Hi all, >> >> I'm experiencing a bizarre issue when trying to implement zero-flux external >> boundary condition when solving for transient solutions. My understanding is >> that it is the default setting in FiPy 3. Indeed, without specifying it, the >> solutions (attached) remains at unity throughout the simulation duration. >> However, when I manually tried to fix the exterior face gradient to zero >> with var.faceGrad.constrain(0, where=m.exteriorFaces), the solution began to >> blow up (as shown in the printed statements). Why does this happen? Is there >> a hidden conflict in conservation settings I should be careful of? >> >> Thanks, >> Yun >> >> >> -- >> Yun Tao >> Postdoc >> Center for Infectious Disease Dynamics >> Pennsylvania State University >> State College, PA 16803 >> <0flux_bc_fipylistserve.py> >> _______________________________________________ >> 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 ] > > > > > -- > Yun Tao > Postdoc > Center for Infectious Disease Dynamics > Pennsylvania State University > State College, PA 16803 > _______________________________________________ > 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 ]
