Sending to the fipy list...

On Wed, Dec 8, 2021 at 10:41 AM Daniel Wheeler
<[email protected]> wrote:
>
> "repr(coeff.faceValue)" is showing you the dependency tree of the
> variable "coeff.faceValue". For example,
>
> ~~~~
> >>> from fipy import Variable
> >>> a = Variable([1, 2])
> >>> b = Variable([2, 3])
> >>> c = a * b
> >>> repr(c)
> '(Variable(value=array([1, 2])) * Variable(value=array([2, 3])))'
> >>> print(c)
> [2 6]
> >>> c.value
> array([2, 6])
> >>> c
> (Variable(value=array([1, 2])) * Variable(value=array([2, 3])))
> ~~~~
>
> So, if we ask for the representation of "c" then we see the variables
> that "c" depends on, not the evaluation of c. Printing "c" directly or
> asking for its value forces an evaluation.
>
> On Wed, Dec 8, 2021 at 8:21 AM Urzainqui Aramburu Iñaki (LUKE)
> <[email protected]> wrote:
> >
> > Hi again,
> >
> > Thanks a lot for your support! This solved my issues.
> >
> > Small comment:
> > After a bit of tinkering, it seems like coeff.faceValue does indeed return 
> > a FaceVariable object. However, if I run print(repr(coeff.faceValue)), it 
> > does not look like it gives a lazy evaluation of the dependent variables. 
> > Instead, I get a bunch of values (numbers) printed out.
> >
> >
> > Iñaki
> > ________________________________
> > From: Daniel Wheeler <[email protected]>
> > Sent: Tuesday, November 30, 2021 09:16 PM
> > To: Urzainqui Aramburu Iñaki (LUKE) <[email protected]>
> > Cc: [email protected] <[email protected]>
> > Subject: Re: [fipy] masked diffusion term, or gmsh physical groups
> >
> > Hi Iñaki,
> >
> > Thanks for your questions. This shouldn't be a problem for FiPy.
> >
> > There are some notable errors in the code you have in the image in
> > your email. Rather than point out the issues I'm just going to provide
> > an example that works.
> >
> > ~~~~
> > from fipy import CellVariable, Grid2D, FaceVariable, DiffusionTerm,
> > Viewer, TransientTerm
> >
> > mesh = Grid2D(nx=2, ny=2, Lx=1.0, Ly=1.0)
> >
> > phi = CellVariable(mesh=mesh, value=0.01, hasOld=True)
> >
> > X, Y = mesh.faceCenters
> >
> > coeff = phi.faceValue * (X != 0.5)
> >
> > print(coeff)
> > # bcs on left side
> > phi.constrain(0.01, where=(X < 0.5) & (Y == 0.0))
> > phi.constrain(1, where=(X < 0.5) & (Y == 1.0))
> >
> > # bcs on right side
> > phi.constrain(1, where=(X > 0.5) & (Y == 0.0))
> > phi.constrain(0.01, where=(X > 0.5) & (Y == 1.0))
> >
> > eqn = (TransientTerm() == DiffusionTerm(coeff))
> >
> > for i in range(10):
> >     phi.updateOld()
> >     eqn.solve(phi, dt=0.1)
> >
> > Viewer(phi).plot()
> >
> > print(phi)
> >
> > input('stop')
> > ~~~~
> >
> > The internal faces in the Y direction have 0 diffusion in the above
> > while the diffusion coefficient is also dependent on phi. The left and
> > right side are therefore independent problems. Here the mask is X !=
> > 0.5. In your example define the diffusion coefficient with
> >
> >     coeff = diff.faceValue * mask
> >
> > diff.faceValue is a face variable. coeff will then depend on diff and
> > diff with depend on phi. You can check that your coeff has the correct
> > dependencies using "print(repr(coeff))" in your code to make
> > absolutely sure you trust the lazy evaluation. That will print the
> > dependency tree. Probably avoid having a function to construct the
> > diffusion coefficient. It shouldn't matter if done correctly, but
> > doesn't buy you anything as it won't be called more than once.
> >
> > I hope that helps.
> >
> > Cheers,
> >
> > Daniel
> >
> > On Thu, Nov 4, 2021 at 5:23 AM Urzainqui Aramburu Iñaki (Luke)
> > <[email protected]> wrote:
> > >
> > > Dear FiPy developers,
> > >
> > > I am trying to solve a diffusion PDE with a non-linear diffusion term in 
> > > a domain that contains 2 subdomains. Those subdomains are not specified 
> > > in the mesh, i.e., the mesh treats everything as a single domain. I want 
> > > to achieve diffusivity = 0 between the faces of one of the subdomains. In 
> > > other words, I am trying to mask the diffusivity FaceVariable so that 
> > > some faces (given, for instance, by the values of a fp.FaceVariable) are 
> > > identically 0.
> > >
> > > My approach so far is the following:
> > >
> > > eq = fp.TransientTerm() == 
> > > fp.DiffusionTerm(coeff=compute_diffusivity(phi))
> > >
> > > And compute_diffusivity is given by:
> > >
> > > def compute_diffusivity(phi):
> > >     diff = fp.numerix.sqrt(2.0*phi) / fp.numerix.exp(4.0*phi**2)
> > >     diff_face = fp.FaceVariable(mesh=mesh, 
> > > value=diff.arithmeticFaceValue.value)
> > >
> > >     # Create mask of subdomain, all ones except zeros where needed
> > >     mask_face = fp.FaceVariable(...)
> > >
> > >     # Apply mask to diffusivity and return.
> > >     # Multiplication works because the target value at the mask is zero.
> > >     return fp.FaceVariable(mesh=mesh, value=diff_face.value * 
> > > mask_face.value)
> > >
> > > FiPy does solve the problem, but I am worried that it does not understand 
> > > that the diffusion variable is non-linear. I suspect that the problem is 
> > > that the .value method used to construct the FaceVariable is determining 
> > > its value directly, without the lazy evaluation of phi.
> > >
> > > How can I solve this problem? What is the best way to set such a mask in 
> > > FiPy?
> > >
> > > Another approach might be to define the subdomains within the mesh (using 
> > > Gmsh 'physical groups', perhaps?) and then reading and using it in FiPy. 
> > > However, I don't know how to do either of those things. (Extra info: I 
> > > need my code to be useful both with meshes created within FiPy and  by 
> > > Gmsh).
> > >
> > > I would really appreciate your help. Thank you in advance!
> > >
> > > Iñaki
> > >
> > > --
> > > To unsubscribe from this group, send email to 
> > > [email protected]
> > >
> > > View this message at https://list.nist.gov/fipy
> > > ---
> > > To unsubscribe from this group and stop receiving emails from it, send an 
> > > email to [email protected].
> >
> >
> >
> > --
> > Daniel Wheeler
>
>
>
> --
> Daniel Wheeler



-- 
Daniel Wheeler

-- 
To unsubscribe from this group, send email to [email protected]

View this message at https://list.nist.gov/fipy
--- 
To unsubscribe from this group and stop receiving emails from it, send an email 
to [email protected].

Reply via email to