Re: Coupling Boundary Condition of one PDE with source term of another PDE
On Tue, Sep 13, 2016 at 12:59 PM, Guyer, Jonathan E. Dr. (Fed) wrote: > I'm not sure; Daniel might have a better idea. > > One thing you can do to manage this (if the anisotropic diffusion is > necessary) is change 2c) to: > >c) eq1 = fp.DiffusionTerm(coeff=S * [[[1., 0], [0., 0.]]]) == f(w, a) > # matrix co-oeff has now been simplified by using a single facevariable whose > values are appropriately zeroed out in the 2D grid Krishna, the above is correct the coefficient can be both a rank 2 cell (or face) variable and have it's value vary in space. You could do my_coeff = fp.CellVariable(mesh=mesh, rank=2, value=0.) and then set the value explicitly my_coeff[0, 0] = x**2 ## or however you want to set the values and then set the coeff in the diffusion term to be "my_coeff". This is exactly equivalent to what Jon is saying and similar to what you suggested. It should work fine. -- 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: Coupling Boundary Condition of one PDE with source term of another PDE
I'm not sure; Daniel might have a better idea. One thing you can do to manage this (if the anisotropic diffusion is necessary) is change 2c) to: c) eq1 = fp.DiffusionTerm(coeff=S * [[[1., 0], [0., 0.]]]) == f(w, a) # matrix co-oeff has now been simplified by using a single facevariable whose values are appropriately zeroed out in the 2D grid > On Sep 12, 2016, at 4:34 PM, Gopalakrishnan, Krishnakumar > wrote: > > Dear all, > > We have a question regarding the original posting to this list (below) > > To recap and quickly summarise, we have > 1. a 1D elliptic PDE which is defined along the x-axis (at the top of a 2D > grid) and > 2. a standard diffusion PDE, which is restricted to diffuse only along the > y-axis. > > The BC of the "vertically diffusing PDE" depends on the values of the axial > PDE solution along the top of the mesh. > > This question is based on Jon's gist code 'crazycoupling2.ipynb' here, > https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d . > > For the axial (elliptic) PDE, to confine the direction of solution, in Input > cell 33 of the jupyter notebook, we can see that the following has been > defined. > > eq1 = fp.DiffusionTerm(coeff=[[[S, 0.], [0., 0.]]]) == f(w, a) > > > Is defining coeff=[[[S, 0.], [0., 0.]]]) the only way to achieve this effect > ? I am asking because, we have a high-order tensor coefficient which is > dependent spatially upon the mesh's x- co-ordinates (as x**2.0), and thus the > value S is not a simple scalar anymore. > > Does the following (simpler) approach work and produce an identical result ? > > 1. declare S as a rank-1 faceVariable in the 2D mesh, initialised to zero in > all faces > 2. Selectively set S only at the vertical faces bordering the top-cells of > the grid using the construct, >a) x = mesh.cellCenters[0] , Y = mesh.faceCenters[1] >b) S.setValue(x**2.0, where = (Y == 1.0 - 0.5*dy)) , where 1.0 is the > height of the mesh, and dy is the constant inter-nodal spacing. >c) eq1 = fp.DiffusionTerm(coeff=S) == f(w, a) # matrix co-oeff has now > been simplified by using a single facevariable whose values are appropriately > zeroed out in the 2D grid > > > In effect, my question is that, if the diffusion coefficient along all other > faces is set to zero, do we still need the coeff=[[[S, 0.], [0., 0.]]]) in > order to restrict the PDE to the axial direction ? > > > Best Regards > > Ian and Krishna > > > > > > -Original Message----- > From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of > Guyer, Jonathan E. Dr. (Fed) > Sent: Thursday, July 14, 2016 5:26 PM > To: FIPY > Subject: Re: Coupling Boundary Condition of one PDE with source term of > another PDE > > My gut reaction is that the second implementation is probably better, as it > should be possible to make things more implicit and coupled, even though > you're "wasting" calculation of PDE1 over most of the 2D domain where you > don't care about it. > >> On Jul 14, 2016, at 11:18 AM, Guyer, Jonathan E. Dr. (Fed) >> wrote: >> >>> What’s the best way to implement this problem in FiPy? >> >> Don't! >> >> >> Assuming you won't take that advise, I've posted a couple of attempts at >> this problem at: >> >> https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d >> >> >> >>> On Jul 9, 2016, at 1:37 PM, Gopalakrishnan, Krishnakumar >>> wrote: >>> >>> Hi, >>> >>> **We’re sorry, the previous posting omitted the complete definition of the >>> boundary condition for the 2nd PDE, and also a couple of typos are fixed >>> below. >>> >>> We have a system of equations, wherein the BC of one PDE couples with the >>> source term of another PDE. >>> >>> We have a regular 2D unit grid in x and y. >>> >>> There are two PDEs to be solved: >>> a)The first PDE (elliptic diffusion problem) is defined only at y = 1, >>> acting along the x-axis (i.e. it acts in the x-direction and only along the >>> top of the cartesian grid). This x-axis is discretized with a fixed >>> grid-spacing, generating a finite number of nodes. Let this set of >>> nodes/co-ordinates be represented by ‘X’. >>> b) The second PDE is a time-varying diffusion problem. This is defined >>> only along the y-axis, but for all x-nodes (i.e. for all ‘X’), where the >>> 1st PDE is being solved. >>> >>>
RE: Coupling Boundary Condition of one PDE with source term of another PDE
Dear all, We have a question regarding the original posting to this list (below) To recap and quickly summarise, we have 1. a 1D elliptic PDE which is defined along the x-axis (at the top of a 2D grid) and 2. a standard diffusion PDE, which is restricted to diffuse only along the y-axis. The BC of the "vertically diffusing PDE" depends on the values of the axial PDE solution along the top of the mesh. This question is based on Jon's gist code 'crazycoupling2.ipynb' here, https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d . For the axial (elliptic) PDE, to confine the direction of solution, in Input cell 33 of the jupyter notebook, we can see that the following has been defined. eq1 = fp.DiffusionTerm(coeff=[[[S, 0.], [0., 0.]]]) == f(w, a) Is defining coeff=[[[S, 0.], [0., 0.]]]) the only way to achieve this effect ? I am asking because, we have a high-order tensor coefficient which is dependent spatially upon the mesh's x- co-ordinates (as x**2.0), and thus the value S is not a simple scalar anymore. Does the following (simpler) approach work and produce an identical result ? 1. declare S as a rank-1 faceVariable in the 2D mesh, initialised to zero in all faces 2. Selectively set S only at the vertical faces bordering the top-cells of the grid using the construct, a) x = mesh.cellCenters[0] , Y = mesh.faceCenters[1] b) S.setValue(x**2.0, where = (Y == 1.0 - 0.5*dy)) , where 1.0 is the height of the mesh, and dy is the constant inter-nodal spacing. c) eq1 = fp.DiffusionTerm(coeff=S) == f(w, a) # matrix co-oeff has now been simplified by using a single facevariable whose values are appropriately zeroed out in the 2D grid In effect, my question is that, if the diffusion coefficient along all other faces is set to zero, do we still need the coeff=[[[S, 0.], [0., 0.]]]) in order to restrict the PDE to the axial direction ? Best Regards Ian and Krishna -Original Message- From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of Guyer, Jonathan E. Dr. (Fed) Sent: Thursday, July 14, 2016 5:26 PM To: FIPY Subject: Re: Coupling Boundary Condition of one PDE with source term of another PDE My gut reaction is that the second implementation is probably better, as it should be possible to make things more implicit and coupled, even though you're "wasting" calculation of PDE1 over most of the 2D domain where you don't care about it. > On Jul 14, 2016, at 11:18 AM, Guyer, Jonathan E. Dr. (Fed) > wrote: > >> What’s the best way to implement this problem in FiPy? > > Don't! > > > Assuming you won't take that advise, I've posted a couple of attempts at this > problem at: > > https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d > > > >> On Jul 9, 2016, at 1:37 PM, Gopalakrishnan, Krishnakumar >> wrote: >> >> Hi, >> >> **We’re sorry, the previous posting omitted the complete definition of the >> boundary condition for the 2nd PDE, and also a couple of typos are fixed >> below. >> >> We have a system of equations, wherein the BC of one PDE couples with the >> source term of another PDE. >> >> We have a regular 2D unit grid in x and y. >> >> There are two PDEs to be solved: >> a)The first PDE (elliptic diffusion problem) is defined only at y = 1, >> acting along the x-axis (i.e. it acts in the x-direction and only along the >> top of the cartesian grid). This x-axis is discretized with a fixed >> grid-spacing, generating a finite number of nodes. Let this set of >> nodes/co-ordinates be represented by ‘X’. >> b) The second PDE is a time-varying diffusion problem. This is defined >> only along the y-axis, but for all x-nodes (i.e. for all ‘X’), where the 1st >> PDE is being solved. >> >> Mathematically (plain-text version) : >> >> PDE1: divergence(S * grad(a)) = f(w,a) >> BC1: a = 0 , at x = 0 (dirichelet) >> BC2: d_a/dx = 1 at x = 1 (Neumann) >> >> D = faceVariable(rank=2, value = 1.0) # Rank 2 tensor for anisotropic >> diffusion (i.e. allowed only along the y-axis) >> >> PDE2: partial_B/partial_t = divergence( [[[0, 0], [0, D]]] * grad(B) ) >> BC1: B = 0, at y = 0 at all ‘X’ (dirichelet), i.e. along the bottom face >> BC2: d_B/dy = g(w,a) (Neumann) at y =1, i.e. along the top face of the grid >> >> f and g are linear functions in ‘w(x,y,t)’ and ‘a(x,y,t)’. B is defined in >> the 2D grid as ‘B(x,y,t)’. >> >> Importantly, w = B at y = 1, for all ‘X’; i.e., the BC2 of the 2nd PDE >> couples with the Implicit Source Term of the 1st PDE along the top face of >> the Cartesian mesh. >> >> What’s the best way to implement this problem in FiPy? >> >> Best Regards, >> >> Krishna ___ fipy mailing list fipy@nist.gov http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
Re: Coupling Boundary Condition of one PDE with source term of another PDE
In your initial condition, w is zero everywhere, correct? In that case, your PDE is indistinguishable from Laplace's equation and so has infinitely many solutions. You might try setting your initial condition of w to something nonzero to see if that helps initial convergence. > On Jul 26, 2016, at 10:10 AM, Gopalakrishnan, Krishnakumar > wrote: > > Hi Jon, > > Thanks a lot for replying. > > Based a reply to one of our earlier posts on this mailing list on 26 May 2016 > *, we understand that: > > If the right-hand-side (i.e. source-term) of the elliptic PDE is a function > of the field-variable, then this uniquely pins down the solution vector > (even though we have only a Neumann BC). I am quoting the relevant sentence > from that post below. > > "However, the form of the source is important here because if the integral of > the source is non-zero, then there is no steady solution. If the integral > _is_ zero and the source is not a function of \phi, then the system admits an > infinite family of solutions, all shifted by a constant." - Raymond Smith > > * http://thread.gmane.org/gmane.comp.python.fipy/4110/focus=4111 > > > > We tried this approach on a standalone elliptic PDE with implicit source term > (albeit without the coupling into the 2D domain with anisotropic diffusion), > and it was successful. However, when we modify your gist code, it throws up > errors with Zero-flux Neumann BC on the left boundary (even though the > f(w,a) is now replaced with ImplicitSourceTerm (coeff=w**3.0). I assume that > this should have uniquely constrained the solution, as it did in our > standalone implementation. Kindly correct us if our understanding is lacking > in this context. > > We definitely need to implement the zero-flux Neumann BC anyway, with the > crazy coupling described in the original post. Can you recommend a solution > for implementing this ? > > > Best Regards, > > Ian & Krishna > > > -Original Message- > From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of > Guyer, Jonathan E. Dr. (Fed) > Sent: Tuesday, July 26, 2016 1:16 PM > To: FIPY > Subject: Re: Coupling Boundary Condition of one PDE with source term of > another PDE > > The source in Poisson's equation determines the curvature of the solution. > Neumann boundary conditions determine the slope at the boundaries. There are > still infinitely many solutions that satisfy those conditions. > >> On Jul 23, 2016, at 7:47 AM, Gopalakrishnan, Krishnakumar >> wrote: >> >> Hi Jon, >> >> Thanks a lot for your valuable code in getting us onto the right track. >> >> However, there is an immediate problem when we start to modify the code >> slightly. If the BC1 of PDE1 is changed from Dirichelet to (zero-flux) >> Neumann, the solver stagnates. In addition, we get the error >> 'ScalarQuantityOutOfRangeWarning: A scalar quantity became too small or too >> large to continue computing. Iterations: 15. Relative error: 1.63296e+16' >> when the interpreter gets to the line, res1 = eq1.sweep(a). The results >> have a lot of NaNs and a few other error and warning messages follow. This >> is the case in both the crazy_coupling1 and crazy_coupling2 codes. >> >> Zero-flux Neumann is a valid BC for PDE1 since the PDE is implicit, i.e. the >> source term on the right is defined by f(w, a), where a is the CellVariable. >> This should pin down the value at the left boundary and hence make the >> physics valid. How can we massage your gist code to be amenable to this? >> >> >> Best Regards, >> >> Krishna >> >> >> -Original Message- >> From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of >> Guyer, Jonathan E. Dr. (Fed) >> Sent: Thursday, July 14, 2016 5:26 PM >> To: FIPY >> Subject: Re: Coupling Boundary Condition of one PDE with source term of >> another PDE >> >> My gut reaction is that the second implementation is probably better, as it >> should be possible to make things more implicit and coupled, even though >> you're "wasting" calculation of PDE1 over most of the 2D domain where you >> don't care about it. >> >>> On Jul 14, 2016, at 11:18 AM, Guyer, Jonathan E. Dr. (Fed) >>> wrote: >>> >>>> What’s the best way to implement this problem in FiPy? >>> >>> Don't! >>> >>> >>> Assuming you won't take that advise, I've posted a couple of attempts at >>> this problem at
RE: Coupling Boundary Condition of one PDE with source term of another PDE
Hi Jon, Thanks a lot for replying. Based a reply to one of our earlier posts on this mailing list on 26 May 2016 *, we understand that: If the right-hand-side (i.e. source-term) of the elliptic PDE is a function of the field-variable, then this uniquely pins down the solution vector (even though we have only a Neumann BC). I am quoting the relevant sentence from that post below. "However, the form of the source is important here because if the integral of the source is non-zero, then there is no steady solution. If the integral _is_ zero and the source is not a function of \phi, then the system admits an infinite family of solutions, all shifted by a constant." - Raymond Smith * http://thread.gmane.org/gmane.comp.python.fipy/4110/focus=4111 We tried this approach on a standalone elliptic PDE with implicit source term (albeit without the coupling into the 2D domain with anisotropic diffusion), and it was successful. However, when we modify your gist code, it throws up errors with Zero-flux Neumann BC on the left boundary (even though the f(w,a) is now replaced with ImplicitSourceTerm (coeff=w**3.0). I assume that this should have uniquely constrained the solution, as it did in our standalone implementation. Kindly correct us if our understanding is lacking in this context. We definitely need to implement the zero-flux Neumann BC anyway, with the crazy coupling described in the original post. Can you recommend a solution for implementing this ? Best Regards, Ian & Krishna -Original Message- From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of Guyer, Jonathan E. Dr. (Fed) Sent: Tuesday, July 26, 2016 1:16 PM To: FIPY Subject: Re: Coupling Boundary Condition of one PDE with source term of another PDE The source in Poisson's equation determines the curvature of the solution. Neumann boundary conditions determine the slope at the boundaries. There are still infinitely many solutions that satisfy those conditions. > On Jul 23, 2016, at 7:47 AM, Gopalakrishnan, Krishnakumar > wrote: > > Hi Jon, > > Thanks a lot for your valuable code in getting us onto the right track. > > However, there is an immediate problem when we start to modify the code > slightly. If the BC1 of PDE1 is changed from Dirichelet to (zero-flux) > Neumann, the solver stagnates. In addition, we get the error > 'ScalarQuantityOutOfRangeWarning: A scalar quantity became too small or too > large to continue computing. Iterations: 15. Relative error: 1.63296e+16' > when the interpreter gets to the line, res1 = eq1.sweep(a). The results have > a lot of NaNs and a few other error and warning messages follow. This is the > case in both the crazy_coupling1 and crazy_coupling2 codes. > > Zero-flux Neumann is a valid BC for PDE1 since the PDE is implicit, i.e. the > source term on the right is defined by f(w, a), where a is the CellVariable. > This should pin down the value at the left boundary and hence make the > physics valid. How can we massage your gist code to be amenable to this? > > > Best Regards, > > Krishna > > > -Original Message- > From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of > Guyer, Jonathan E. Dr. (Fed) > Sent: Thursday, July 14, 2016 5:26 PM > To: FIPY > Subject: Re: Coupling Boundary Condition of one PDE with source term of > another PDE > > My gut reaction is that the second implementation is probably better, as it > should be possible to make things more implicit and coupled, even though > you're "wasting" calculation of PDE1 over most of the 2D domain where you > don't care about it. > >> On Jul 14, 2016, at 11:18 AM, Guyer, Jonathan E. Dr. (Fed) >> wrote: >> >>> What’s the best way to implement this problem in FiPy? >> >> Don't! >> >> >> Assuming you won't take that advise, I've posted a couple of attempts at >> this problem at: >> >> https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d >> >> >> >>> On Jul 9, 2016, at 1:37 PM, Gopalakrishnan, Krishnakumar >>> wrote: >>> >>> Hi, >>> >>> **We’re sorry, the previous posting omitted the complete definition of the >>> boundary condition for the 2nd PDE, and also a couple of typos are fixed >>> below. >>> >>> We have a system of equations, wherein the BC of one PDE couples with the >>> source term of another PDE. >>> >>> We have a regular 2D unit grid in x and y. >>> >>> There are two PDEs to be solved: >>> a)The first PDE (elliptic diffusion problem) is defined only at y = 1, >>> acting a
Re: Coupling Boundary Condition of one PDE with source term of another PDE
The source in Poisson's equation determines the curvature of the solution. Neumann boundary conditions determine the slope at the boundaries. There are still infinitely many solutions that satisfy those conditions. > On Jul 23, 2016, at 7:47 AM, Gopalakrishnan, Krishnakumar > wrote: > > Hi Jon, > > Thanks a lot for your valuable code in getting us onto the right track. > > However, there is an immediate problem when we start to modify the code > slightly. If the BC1 of PDE1 is changed from Dirichelet to (zero-flux) > Neumann, the solver stagnates. In addition, we get the error > 'ScalarQuantityOutOfRangeWarning: A scalar quantity became too small or too > large to continue computing. Iterations: 15. Relative error: 1.63296e+16' > when the interpreter gets to the line, res1 = eq1.sweep(a). The results have > a lot of NaNs and a few other error and warning messages follow. This is the > case in both the crazy_coupling1 and crazy_coupling2 codes. > > Zero-flux Neumann is a valid BC for PDE1 since the PDE is implicit, i.e. the > source term on the right is defined by f(w, a), where a is the CellVariable. > This should pin down the value at the left boundary and hence make the > physics valid. How can we massage your gist code to be amenable to this? > > > Best Regards, > > Krishna > > > -Original Message- > From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of > Guyer, Jonathan E. Dr. (Fed) > Sent: Thursday, July 14, 2016 5:26 PM > To: FIPY > Subject: Re: Coupling Boundary Condition of one PDE with source term of > another PDE > > My gut reaction is that the second implementation is probably better, as it > should be possible to make things more implicit and coupled, even though > you're "wasting" calculation of PDE1 over most of the 2D domain where you > don't care about it. > >> On Jul 14, 2016, at 11:18 AM, Guyer, Jonathan E. Dr. (Fed) >> wrote: >> >>> What’s the best way to implement this problem in FiPy? >> >> Don't! >> >> >> Assuming you won't take that advise, I've posted a couple of attempts at >> this problem at: >> >> https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d >> >> >> >>> On Jul 9, 2016, at 1:37 PM, Gopalakrishnan, Krishnakumar >>> wrote: >>> >>> Hi, >>> >>> **We’re sorry, the previous posting omitted the complete definition of the >>> boundary condition for the 2nd PDE, and also a couple of typos are fixed >>> below. >>> >>> We have a system of equations, wherein the BC of one PDE couples with the >>> source term of another PDE. >>> >>> We have a regular 2D unit grid in x and y. >>> >>> There are two PDEs to be solved: >>> a)The first PDE (elliptic diffusion problem) is defined only at y = 1, >>> acting along the x-axis (i.e. it acts in the x-direction and only along the >>> top of the cartesian grid). This x-axis is discretized with a fixed >>> grid-spacing, generating a finite number of nodes. Let this set of >>> nodes/co-ordinates be represented by ‘X’. >>> b) The second PDE is a time-varying diffusion problem. This is defined >>> only along the y-axis, but for all x-nodes (i.e. for all ‘X’), where the >>> 1st PDE is being solved. >>> >>> Mathematically (plain-text version) : >>> >>> PDE1: divergence(S * grad(a)) = f(w,a) >>> BC1: a = 0 , at x = 0 (dirichelet) >>> BC2: d_a/dx = 1 at x = 1 (Neumann) >>> >>> D = faceVariable(rank=2, value = 1.0) # Rank 2 tensor for anisotropic >>> diffusion (i.e. allowed only along the y-axis) >>> >>> PDE2: partial_B/partial_t = divergence( [[[0, 0], [0, D]]] * grad(B) ) >>> BC1: B = 0, at y = 0 at all ‘X’ (dirichelet), i.e. along the bottom face >>> BC2: d_B/dy = g(w,a) (Neumann) at y =1, i.e. along the top face of the grid >>> >>> f and g are linear functions in ‘w(x,y,t)’ and ‘a(x,y,t)’. B is defined in >>> the 2D grid as ‘B(x,y,t)’. >>> >>> Importantly, w = B at y = 1, for all ‘X’; i.e., the BC2 of the 2nd PDE >>> couples with the Implicit Source Term of the 1st PDE along the top face of >>> the Cartesian mesh. >>> >>> What’s the best way to implement this problem in FiPy? >>> >>> Best Regards, >>> >>> Krishna >>> >>> >>> PS: Since
RE: Coupling Boundary Condition of one PDE with source term of another PDE
Hi Jon, Thanks a lot for your valuable code in getting us onto the right track. However, there is an immediate problem when we start to modify the code slightly. If the BC1 of PDE1 is changed from Dirichelet to (zero-flux) Neumann, the solver stagnates. In addition, we get the error 'ScalarQuantityOutOfRangeWarning: A scalar quantity became too small or too large to continue computing. Iterations: 15. Relative error: 1.63296e+16' when the interpreter gets to the line, res1 = eq1.sweep(a). The results have a lot of NaNs and a few other error and warning messages follow. This is the case in both the crazy_coupling1 and crazy_coupling2 codes. Zero-flux Neumann is a valid BC for PDE1 since the PDE is implicit, i.e. the source term on the right is defined by f(w, a), where a is the CellVariable. This should pin down the value at the left boundary and hence make the physics valid. How can we massage your gist code to be amenable to this? Best Regards, Krishna -Original Message- From: fipy-boun...@nist.gov [mailto:fipy-boun...@nist.gov] On Behalf Of Guyer, Jonathan E. Dr. (Fed) Sent: Thursday, July 14, 2016 5:26 PM To: FIPY Subject: Re: Coupling Boundary Condition of one PDE with source term of another PDE My gut reaction is that the second implementation is probably better, as it should be possible to make things more implicit and coupled, even though you're "wasting" calculation of PDE1 over most of the 2D domain where you don't care about it. > On Jul 14, 2016, at 11:18 AM, Guyer, Jonathan E. Dr. (Fed) > wrote: > >> What’s the best way to implement this problem in FiPy? > > Don't! > > > Assuming you won't take that advise, I've posted a couple of attempts at this > problem at: > > https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d > > > >> On Jul 9, 2016, at 1:37 PM, Gopalakrishnan, Krishnakumar >> wrote: >> >> Hi, >> >> **We’re sorry, the previous posting omitted the complete definition of the >> boundary condition for the 2nd PDE, and also a couple of typos are fixed >> below. >> >> We have a system of equations, wherein the BC of one PDE couples with the >> source term of another PDE. >> >> We have a regular 2D unit grid in x and y. >> >> There are two PDEs to be solved: >> a)The first PDE (elliptic diffusion problem) is defined only at y = 1, >> acting along the x-axis (i.e. it acts in the x-direction and only along the >> top of the cartesian grid). This x-axis is discretized with a fixed >> grid-spacing, generating a finite number of nodes. Let this set of >> nodes/co-ordinates be represented by ‘X’. >> b) The second PDE is a time-varying diffusion problem. This is defined >> only along the y-axis, but for all x-nodes (i.e. for all ‘X’), where the 1st >> PDE is being solved. >> >> Mathematically (plain-text version) : >> >> PDE1: divergence(S * grad(a)) = f(w,a) >> BC1: a = 0 , at x = 0 (dirichelet) >> BC2: d_a/dx = 1 at x = 1 (Neumann) >> >> D = faceVariable(rank=2, value = 1.0) # Rank 2 tensor for anisotropic >> diffusion (i.e. allowed only along the y-axis) >> >> PDE2: partial_B/partial_t = divergence( [[[0, 0], [0, D]]] * grad(B) ) >> BC1: B = 0, at y = 0 at all ‘X’ (dirichelet), i.e. along the bottom face >> BC2: d_B/dy = g(w,a) (Neumann) at y =1, i.e. along the top face of the grid >> >> f and g are linear functions in ‘w(x,y,t)’ and ‘a(x,y,t)’. B is defined in >> the 2D grid as ‘B(x,y,t)’. >> >> Importantly, w = B at y = 1, for all ‘X’; i.e., the BC2 of the 2nd PDE >> couples with the Implicit Source Term of the 1st PDE along the top face of >> the Cartesian mesh. >> >> What’s the best way to implement this problem in FiPy? >> >> Best Regards, >> >> Krishna >> >> >> PS: Since the problem is mathematically harder to express in plain-text, >> here is a HTML formatted version of the PDEs and BCs. >> PDE1: >> BC1:(Dirichelet) >> BC2:(Neumann) >> >> D = faceVariable(mesh=pos_p2d_mesh, rank=2, value = 1.0) # Rank 2 >> tensor for anisotropic diffusion (i.e. allowed only along the y-axis) >> >> PDE2: = >> BC1 is: i.e. along the bottom face of the grid >> BC2 is: (Neumann) at y =1, i.e. along the top face of the >> grid >> >> f and g are linear functions in and. >> is defined in the 2D grid as >> >> Importantly, i.e., the BC2 of the 2nd PDE couples with the >> Implicit Source Term of the 1st PDE along the top face of the Cartesian mesh. >
Re: Coupling Boundary Condition of one PDE with source term of another PDE
My gut reaction is that the second implementation is probably better, as it should be possible to make things more implicit and coupled, even though you're "wasting" calculation of PDE1 over most of the 2D domain where you don't care about it. > On Jul 14, 2016, at 11:18 AM, Guyer, Jonathan E. Dr. (Fed) > wrote: > >> What’s the best way to implement this problem in FiPy? > > Don't! > > > Assuming you won't take that advise, I've posted a couple of attempts at this > problem at: > > https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d > > > >> On Jul 9, 2016, at 1:37 PM, Gopalakrishnan, Krishnakumar >> wrote: >> >> Hi, >> >> **We’re sorry, the previous posting omitted the complete definition of the >> boundary condition for the 2nd PDE, and also a couple of typos are fixed >> below. >> >> We have a system of equations, wherein the BC of one PDE couples with the >> source term of another PDE. >> >> We have a regular 2D unit grid in x and y. >> >> There are two PDEs to be solved: >> a)The first PDE (elliptic diffusion problem) is defined only at y = 1, >> acting along the x-axis (i.e. it acts in the x-direction and only along the >> top of the cartesian grid). This x-axis is discretized with a fixed >> grid-spacing, generating a finite number of nodes. Let this set of >> nodes/co-ordinates be represented by ‘X’. >> b) The second PDE is a time-varying diffusion problem. This is defined >> only along the y-axis, but for all x-nodes (i.e. for all ‘X’), where the 1st >> PDE is being solved. >> >> Mathematically (plain-text version) : >> >> PDE1: divergence(S * grad(a)) = f(w,a) >> BC1: a = 0 , at x = 0 (dirichelet) >> BC2: d_a/dx = 1 at x = 1 (Neumann) >> >> D = faceVariable(rank=2, value = 1.0) # Rank 2 tensor for anisotropic >> diffusion (i.e. allowed only along the y-axis) >> >> PDE2: partial_B/partial_t = divergence( [[[0, 0], [0, D]]] * grad(B) ) >> BC1: B = 0, at y = 0 at all ‘X’ (dirichelet), i.e. along the bottom face >> BC2: d_B/dy = g(w,a) (Neumann) at y =1, i.e. along the top face of the grid >> >> f and g are linear functions in ‘w(x,y,t)’ and ‘a(x,y,t)’. B is defined in >> the 2D grid as ‘B(x,y,t)’. >> >> Importantly, w = B at y = 1, for all ‘X’; i.e., the BC2 of the 2nd PDE >> couples with the Implicit Source Term of the 1st PDE along the top face of >> the Cartesian mesh. >> >> What’s the best way to implement this problem in FiPy? >> >> Best Regards, >> >> Krishna >> >> >> PS: Since the problem is mathematically harder to express in plain-text, >> here is a HTML formatted version of the PDEs and BCs. >> PDE1: >> BC1:(Dirichelet) >> BC2:(Neumann) >> >> D = faceVariable(mesh=pos_p2d_mesh, rank=2, value = 1.0) # Rank 2 >> tensor for anisotropic diffusion (i.e. allowed only along the y-axis) >> >> PDE2: = >> BC1 is: i.e. along the bottom face of the grid >> BC2 is: (Neumann) at y =1, i.e. along the top face of the >> grid >> >> f and g are linear functions in and. >> is defined in the 2D grid as >> >> Importantly, i.e., the BC2 of the 2nd PDE couples with the >> Implicit Source Term of the 1st PDE along the top face of the Cartesian mesh. >> >> ___ >> 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: Coupling Boundary Condition of one PDE with source term of another PDE
> What’s the best way to implement this problem in FiPy? Don't! Assuming you won't take that advise, I've posted a couple of attempts at this problem at: https://gist.github.com/guyer/bb199559c00f6047d466daa18554d83d > On Jul 9, 2016, at 1:37 PM, Gopalakrishnan, Krishnakumar > wrote: > > Hi, > > **We’re sorry, the previous posting omitted the complete definition of the > boundary condition for the 2nd PDE, and also a couple of typos are fixed > below. > > We have a system of equations, wherein the BC of one PDE couples with the > source term of another PDE. > > We have a regular 2D unit grid in x and y. > > There are two PDEs to be solved: > a)The first PDE (elliptic diffusion problem) is defined only at y = 1, > acting along the x-axis (i.e. it acts in the x-direction and only along the > top of the cartesian grid). This x-axis is discretized with a fixed > grid-spacing, generating a finite number of nodes. Let this set of > nodes/co-ordinates be represented by ‘X’. > b) The second PDE is a time-varying diffusion problem. This is defined only > along the y-axis, but for all x-nodes (i.e. for all ‘X’), where the 1st PDE > is being solved. > > Mathematically (plain-text version) : > > PDE1: divergence(S * grad(a)) = f(w,a) > BC1: a = 0 , at x = 0 (dirichelet) > BC2: d_a/dx = 1 at x = 1 (Neumann) > > D = faceVariable(rank=2, value = 1.0) # Rank 2 tensor for anisotropic > diffusion (i.e. allowed only along the y-axis) > > PDE2: partial_B/partial_t = divergence( [[[0, 0], [0, D]]] * grad(B) ) > BC1: B = 0, at y = 0 at all ‘X’ (dirichelet), i.e. along the bottom face > BC2: d_B/dy = g(w,a) (Neumann) at y =1, i.e. along the top face of the grid > > f and g are linear functions in ‘w(x,y,t)’ and ‘a(x,y,t)’. B is defined in > the 2D grid as ‘B(x,y,t)’. > > Importantly, w = B at y = 1, for all ‘X’; i.e., the BC2 of the 2nd PDE > couples with the Implicit Source Term of the 1st PDE along the top face of > the Cartesian mesh. > > What’s the best way to implement this problem in FiPy? > > Best Regards, > > Krishna > > > PS: Since the problem is mathematically harder to express in plain-text, here > is a HTML formatted version of the PDEs and BCs. > PDE1: > BC1:(Dirichelet) > BC2:(Neumann) > > D = faceVariable(mesh=pos_p2d_mesh, rank=2, value = 1.0) # Rank 2 > tensor for anisotropic diffusion (i.e. allowed only along the y-axis) > > PDE2: = > BC1 is: i.e. along the bottom face of the grid > BC2 is: (Neumann) at y =1, i.e. along the top face of the grid > > f and g are linear functions in and. > is defined in the 2D grid as > > Importantly, i.e., the BC2 of the 2nd PDE couples with the > Implicit Source Term of the 1st PDE along the top face of the Cartesian mesh. > > ___ > 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 ]
