Re: Setting up 3D Diffusion

2019-04-04 Thread Guyer, Jonathan E. Dr. (Fed) via fipy
> On Apr 3, 2019, at 6:36 PM, Dario Panada  wrote:
> 
> Many thanks, I think that is making more sense now. I assume I could also 
> change D to be expressed in terms of 0.042**(2./3) to keep consistency? (It 
> would be a bit of a bizarre setup I admit..)

If you choose to scale things on a length of 0.042**(2./3), then you should do 
so for both d(x|y|z) and for D. 

I reiterate, there is no computational benefit to doing so. Further, if it 
matters to you that your cells are 0.042 mm**3, it is unlikely that you'll 
derive any conceptual benefit from a dimensionless representation.

> As a separate question, if not boundary conditions are specified what does 
> this default to? (If anything)

No flux


___
fipy mailing list
fipy@nist.gov
http://www.ctcms.nist.gov/fipy
  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]


Re: Setting up 3D Diffusion

2019-04-03 Thread Dario Panada
Many thanks, I think that is making more sense now. I assume I could also
change D to be expressed in terms of 0.042**(2./3) to keep consistency? (It
would be a bit of a bizarre setup I admit..)

As a separate question, if not boundary conditions are specified what does
this default to? (If anything)


On Wed, Apr 3, 2019 at 10:52 PM Guyer, Jonathan E. Dr. (Fed) via fipy <
fipy@nist.gov> wrote:

> Dario -
>
> 1) FiPy does have a mechanism for specifying units. See:
>
>
> https://www.ctcms.nist.gov/fipy/fipy/generated/fipy.variables.html#fipy.variables.variable.Variable
>
> https://www.ctcms.nist.gov/fipy/fipy/generated/fipy.variables.html#fipy.variables.variable.Variable.inBaseUnits
>
> https://www.ctcms.nist.gov/fipy/fipy/generated/fipy.variables.html#fipy.variables.variable.Variable.inUnitsOf
>
> However, they don't propagate through the solvers, so you end up needing
> to express your equations with dimensionless coefficients composed of
> dimensional Variables. I'm not sure that anybody but me has ever tried to
> use this feature, and I don't use it.
>
> I was not suggesting truly dimensional coefficients, though. Just pick a
> self-consistent system, e.g.,
>
> dx = dy = dz = .042**(1./3) # mm
> dt = 15.125 # s
> D = 0.001 # mm2/s
>
> FiPy doesn't care. You're going to pay the price for multiplying and
> dividing everything by 1, anyway, so you might as well work with terms that
> are easier to think about, at least for me.
>
> 2) I did not mean to suggest that the time step doesn't matter. Only that
> the Fourier limit doesn't apply.
> The Fourier limit is associated with *stability*. An explicit
> discretization "blows up" if the Fourier limit is exceeded, but an implicit
> discretization does not. Accuracy may suffer, however.
>
> This is all covered pretty methodically in examples/diffusion/mesh1D.py.
>
> > On Apr 3, 2019, at 4:35 PM, Dario Panada  wrote:
> >
> > Hi Jonathan,
> >
> > Many thanks for your quick response.
> >
> > Firstly apologies if I'm asking rather basic questions, I'm still new to
> FiPy and learning. I've read the manual but if you think there are any
> sections I would benefit from re-reading, or any additional resources,
> please feel free to point me to them as an answer.
> >
> > 1) So how do I specify/enforce units? For example, I'm setting my
> diffusivity in mm2/s. I could convert this to cm2/s and put that obviously
> different value for D. How does FiPy know I've changed the units? Or should
> I make my units consistent with the spatial discretization? (Ie: If I use
> mm3 then my diffusivity should use mm2.) But then how about the time unit?
> What if I converted from mm2/s to mm2/hour?
> >
> > 2) In my case, changing dt does affect the final output all other things
> being equal. (However, I am using integer for dx etc. Maybe that's why?)
> >
> > Kind Regards,
> > Dario
> >
> > On Wed, Apr 3, 2019 at 8:48 PM Guyer, Jonathan E. Dr. (Fed) via fipy <
> fipy@nist.gov> wrote:
> > Dario -
> >
> > - FiPy's DiffusionTerm is implicit, so there is no Fourier limit on your
> time step.
> >
> > - There's no particular advantage to scaling length and time to 1
> >
> > - Be sure to set dx = dy = dz = 1., a float, not dx = dy = dz = 1, an
> integer
> >
> > - Jon
> >
> > > On Apr 3, 2019, at 11:28 AM, Dario Panada 
> wrote:
> > >
> > > Hello,
> > >
> > > I am trying to use FiPy to build a system to simulate oxygen diffusion
> in tissues in 3D.
> > >
> > > My grid is a 20x20x20 mesh, defined as:
> > >
> > > mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> > >
> > > Where nx = ny = nz = 20 and dx = dy = dz = 1. This approach is taken
> from examples.diffusion.mesh1D, where nx was set to 50 and dx to 1 to
> represent unity.
> > >
> > > I setup my equation as:
> > >
> > > eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> > >
> > > My oxygen diffusion coefficient (D) is 0.001 mm2/s. Each cell is my
> grid corresponds to 0.042 mm3. So (Δx)^3 = 0.042mm^3 which implies (Δx^2)
> ~=0.121 mm2.
> > >
> > > I understand that for 3D problems the Fourier number (Fo) has to be
> less than or equal to 1/8.
> > >
> > > With Fo = DΔt/(Δx^2)  <= 1/8 so Δt <= (Δx^2)/(8D) where replacing
> values I get  Δt <= 15.125 .
> > >
> > > From the above, would it be sensible to setup my script as?
> > > ---
> > > nx = ny = nz = 20
> > > dx = dy = dz = 1
> > > dt = 1
> > > D = 0.121
> > > mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> > >
> > > phi=CellVariable(name="solutionvariable",mesh=mesh,value=0.)
> > >
> > > eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> > >
> > > eq.solve(var=phi,dt=dt)
> > > ---
> > >
> > > Or am I missing something/have interpreted something wrong?
> > >
> > > Kind Regards,
> > > Dario
> > > ___
> > > fipy mailing list
> > > fipy@nist.gov
> > > http://www.ctcms.nist.gov/fipy
> > >  [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
> >
> >
> > 

Re: Setting up 3D Diffusion

2019-04-03 Thread Guyer, Jonathan E. Dr. (Fed) via fipy
Dario -

1) FiPy does have a mechanism for specifying units. See:

https://www.ctcms.nist.gov/fipy/fipy/generated/fipy.variables.html#fipy.variables.variable.Variable
https://www.ctcms.nist.gov/fipy/fipy/generated/fipy.variables.html#fipy.variables.variable.Variable.inBaseUnits
https://www.ctcms.nist.gov/fipy/fipy/generated/fipy.variables.html#fipy.variables.variable.Variable.inUnitsOf

However, they don't propagate through the solvers, so you end up needing to 
express your equations with dimensionless coefficients composed of dimensional 
Variables. I'm not sure that anybody but me has ever tried to use this feature, 
and I don't use it.

I was not suggesting truly dimensional coefficients, though. Just pick a 
self-consistent system, e.g.,

dx = dy = dz = .042**(1./3) # mm
dt = 15.125 # s
D = 0.001 # mm2/s

FiPy doesn't care. You're going to pay the price for multiplying and dividing 
everything by 1, anyway, so you might as well work with terms that are easier 
to think about, at least for me.

2) I did not mean to suggest that the time step doesn't matter. Only that the 
Fourier limit doesn't apply. 
The Fourier limit is associated with *stability*. An explicit discretization 
"blows up" if the Fourier limit is exceeded, but an implicit discretization 
does not. Accuracy may suffer, however. 

This is all covered pretty methodically in examples/diffusion/mesh1D.py.

> On Apr 3, 2019, at 4:35 PM, Dario Panada  wrote:
> 
> Hi Jonathan,
> 
> Many thanks for your quick response.
> 
> Firstly apologies if I'm asking rather basic questions, I'm still new to FiPy 
> and learning. I've read the manual but if you think there are any sections I 
> would benefit from re-reading, or any additional resources, please feel free 
> to point me to them as an answer.
> 
> 1) So how do I specify/enforce units? For example, I'm setting my diffusivity 
> in mm2/s. I could convert this to cm2/s and put that obviously different 
> value for D. How does FiPy know I've changed the units? Or should I make my 
> units consistent with the spatial discretization? (Ie: If I use mm3 then my 
> diffusivity should use mm2.) But then how about the time unit? What if I 
> converted from mm2/s to mm2/hour?
> 
> 2) In my case, changing dt does affect the final output all other things 
> being equal. (However, I am using integer for dx etc. Maybe that's why?)
> 
> Kind Regards,
> Dario
> 
> On Wed, Apr 3, 2019 at 8:48 PM Guyer, Jonathan E. Dr. (Fed) via fipy 
>  wrote:
> Dario -
> 
> - FiPy's DiffusionTerm is implicit, so there is no Fourier limit on your time 
> step.
> 
> - There's no particular advantage to scaling length and time to 1
> 
> - Be sure to set dx = dy = dz = 1., a float, not dx = dy = dz = 1, an integer
> 
> - Jon
> 
> > On Apr 3, 2019, at 11:28 AM, Dario Panada  wrote:
> > 
> > Hello,
> > 
> > I am trying to use FiPy to build a system to simulate oxygen diffusion in 
> > tissues in 3D.
> > 
> > My grid is a 20x20x20 mesh, defined as:
> > 
> > mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> > 
> > Where nx = ny = nz = 20 and dx = dy = dz = 1. This approach is taken from 
> > examples.diffusion.mesh1D, where nx was set to 50 and dx to 1 to represent 
> > unity.
> > 
> > I setup my equation as:
> > 
> > eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> > 
> > My oxygen diffusion coefficient (D) is 0.001 mm2/s. Each cell is my grid 
> > corresponds to 0.042 mm3. So (Δx)^3 = 0.042mm^3 which implies (Δx^2) 
> > ~=0.121 mm2.
> > 
> > I understand that for 3D problems the Fourier number (Fo) has to be less 
> > than or equal to 1/8.
> > 
> > With Fo = DΔt/(Δx^2)  <= 1/8 so Δt <= (Δx^2)/(8D) where replacing values I 
> > get  Δt <= 15.125 .
> > 
> > From the above, would it be sensible to setup my script as?
> > ---
> > nx = ny = nz = 20
> > dx = dy = dz = 1
> > dt = 1
> > D = 0.121
> > mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> > 
> > phi=CellVariable(name="solutionvariable",mesh=mesh,value=0.)
> > 
> > eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> > 
> > eq.solve(var=phi,dt=dt)
> > ---
> > 
> > Or am I missing something/have interpreted something wrong?
> > 
> > Kind Regards,
> > Dario
> > ___
> > 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: Setting up 3D Diffusion

2019-04-03 Thread Dario Panada
Hi Jonathan,

Many thanks for your quick response.

Firstly apologies if I'm asking rather basic questions, I'm still new to
FiPy and learning. I've read the manual but if you think there are any
sections I would benefit from re-reading, or any additional resources,
please feel free to point me to them as an answer.

1) So how do I specify/enforce units? For example, I'm setting my
diffusivity in mm2/s. I could convert this to cm2/s and put that obviously
different value for D. How does FiPy know I've changed the units? Or should
I make my units consistent with the spatial discretization? (Ie: If I use
mm3 then my diffusivity should use mm2.) But then how about the time unit?
What if I converted from mm2/s to mm2/hour?

2) In my case, changing dt *does* affect the final output all other things
being equal. (However, I am using integer for dx etc. Maybe that's why?)

Kind Regards,
Dario

On Wed, Apr 3, 2019 at 8:48 PM Guyer, Jonathan E. Dr. (Fed) via fipy <
fipy@nist.gov> wrote:

> Dario -
>
> - FiPy's DiffusionTerm is implicit, so there is no Fourier limit on your
> time step.
>
> - There's no particular advantage to scaling length and time to 1
>
> - Be sure to set dx = dy = dz = 1., a float, not dx = dy = dz = 1, an
> integer
>
> - Jon
>
> > On Apr 3, 2019, at 11:28 AM, Dario Panada 
> wrote:
> >
> > Hello,
> >
> > I am trying to use FiPy to build a system to simulate oxygen diffusion
> in tissues in 3D.
> >
> > My grid is a 20x20x20 mesh, defined as:
> >
> > mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> >
> > Where nx = ny = nz = 20 and dx = dy = dz = 1. This approach is taken
> from examples.diffusion.mesh1D, where nx was set to 50 and dx to 1 to
> represent unity.
> >
> > I setup my equation as:
> >
> > eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> >
> > My oxygen diffusion coefficient (D) is 0.001 mm2/s. Each cell is my grid
> corresponds to 0.042 mm3. So (Δx)^3 = 0.042mm^3 which implies (Δx^2)
> ~=0.121 mm2.
> >
> > I understand that for 3D problems the Fourier number (Fo) has to be less
> than or equal to 1/8.
> >
> > With Fo = DΔt/(Δx^2)  <= 1/8 so Δt <= (Δx^2)/(8D) where replacing values
> I get  Δt <= 15.125 .
> >
> > From the above, would it be sensible to setup my script as?
> > ---
> > nx = ny = nz = 20
> > dx = dy = dz = 1
> > dt = 1
> > D = 0.121
> > mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> >
> > phi=CellVariable(name="solutionvariable",mesh=mesh,value=0.)
> >
> > eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> >
> > eq.solve(var=phi,dt=dt)
> > ---
> >
> > Or am I missing something/have interpreted something wrong?
> >
> > Kind Regards,
> > Dario
> > ___
> > 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: Setting up 3D Diffusion

2019-04-03 Thread Guyer, Jonathan E. Dr. (Fed) via fipy
Dario -

- FiPy's DiffusionTerm is implicit, so there is no Fourier limit on your time 
step.

- There's no particular advantage to scaling length and time to 1

- Be sure to set dx = dy = dz = 1., a float, not dx = dy = dz = 1, an integer

- Jon

> On Apr 3, 2019, at 11:28 AM, Dario Panada  wrote:
> 
> Hello,
> 
> I am trying to use FiPy to build a system to simulate oxygen diffusion in 
> tissues in 3D.
> 
> My grid is a 20x20x20 mesh, defined as:
> 
> mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> 
> Where nx = ny = nz = 20 and dx = dy = dz = 1. This approach is taken from 
> examples.diffusion.mesh1D, where nx was set to 50 and dx to 1 to represent 
> unity.
> 
> I setup my equation as:
> 
> eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> 
> My oxygen diffusion coefficient (D) is 0.001 mm2/s. Each cell is my grid 
> corresponds to 0.042 mm3. So (Δx)^3 = 0.042mm^3 which implies (Δx^2) ~=0.121 
> mm2.
> 
> I understand that for 3D problems the Fourier number (Fo) has to be less than 
> or equal to 1/8.
> 
> With Fo = DΔt/(Δx^2)  <= 1/8 so Δt <= (Δx^2)/(8D) where replacing values I 
> get  Δt <= 15.125 .
> 
> From the above, would it be sensible to setup my script as?
> ---
> nx = ny = nz = 20
> dx = dy = dz = 1
> dt = 1
> D = 0.121
> mesh=Grid3D(dx=dx,dy=dy,nx=nx,ny=ny, dz=dz, nz=nz)
> 
> phi=CellVariable(name="solutionvariable",mesh=mesh,value=0.)
> 
> eq=TransientTerm()==DiffusionTerm(coeff=D) + sourceGrid - sinkGrid
> 
> eq.solve(var=phi,dt=dt)
> ---
> 
> Or am I missing something/have interpreted something wrong?
> 
> Kind Regards,
> Dario
> ___
> 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 ]