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 <FIPY@nist.gov> 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) > <jonathan.gu...@nist.gov> 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 >> <k.gopalakrishna...@imperial.ac.uk> 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 ]
