Re: Question regarding Boundary Condition Implementation
Jon, Thank you again for your advice. I'm guessing you discovered the mesh problem by visualizing the grid. I will have to remember that for future projects. I will also strongly consider a new boundary condition for the center line of the cylinder. Thank you, Dan On Mon, Apr 22, 2019 at 11:57 AM Guyer, Jonathan E. Dr. (Fed) via fipy < fipy@nist.gov> wrote: > Daniel - > > Two things going on here: > > - Your constraint is being applied at the tip of the wedge of the > cylindrical coordinates. Since this "face" has no area, that Dirichlet > condition doesn't generate any flux in or out of the domain. Offsetting the > cylinder slightly from the origin resolves this: > > mesh = CylindricalGrid1D(Lr = L,nr=nx) + [[L/nx]] > > I recommend giving some thought as to whether it even makes sense to apply > a Dirichlet condition along the centerline of a cylinder. > > - Because you've introduced coefficients which change with temperature, > your equation is now nonlinear and so sweeping is required (2 is apparently > enough). At each time step, execute: > > for sweep in range(2): > res = eq.sweep(dt=dt) > > rho.setValue(rho_l,where = PCM & (T>350) ) > k.setValue(k_l,where = PCM & (T>350)) > Cp.setValue(Cp_l,where = PCM & (T>350)) > > With these changes, the evolution seems generally similar to your Grid1D > case. > > - Jon > > > On Apr 15, 2019, at 1:28 PM, Daniel DeSantis wrote: > > > > Jonathan, > > > > I have just a few follow up questions. Namely, regarding to working > within cylindrical coordinates. If I switch the mesh to cylindrical, it > seems the temperature profile never changes, unless I go to extremely small > time steps. That strikes me as odd. What do you think? > > > > Also, the shifting profile seems to reverse, with the left boundary > condition dropping to ~300K over time while in Cartesian coordinates, the > temperature rises to 600K over time (a result I would expect). I've > included a copy of the work with the changes you've suggested previously. > If you look at lines 23/24, you can toggle the cylindrical grid on and off, > and in lines 30-34, I've included options for the time steps that work for > each respective coordinate system. Do you have any suggestions? > > > > Finally, in reference to the velocity I would like to incorporate, I > understand the difference between Cell Variables and Face Variables, but > why is the rank of the velocity -1? What does the Face Variable rank > indicate? > > > > Thank you, > > Dan > > > > On Thu, Apr 11, 2019 at 8:14 AM Guyer, Jonathan E. Dr. (Fed) via fipy < > fipy@nist.gov> wrote: > > > > > > > On Apr 10, 2019, at 4:42 PM, Daniel DeSantis > wrote: > > > > > > 1) Set a different Cp, rho, and k in just the PCM domain when the > temperature in that domain crosses over a certain melt temperature (350K). > I tried an if statement inside the solving loop and I think that has > worked, based on some print testing I did in the loop. I just wanted to get > your opinion and see if that was the right way to go. > > > > Seems reasonable, as long as you're calling `.setValue()` and not > redefining Cp, rho, and k. > > > > In fact, I wouldn't use an `if`: > > > > >>> for step in range(steps): > > ... T.updateOld() > > ... eq.solve(dt=dt) > > ... rho.setValue(rho_melt, where=PCM & (T > 350)) > > > > > > > 2) I'd like to move this to cylindrical coordinates with the lengths > being radii. I'm not sure if the cylindrical grid is working though. And if > it isn't, is there a way to convert the diffusion equation into something > that would solve the problem in Cartesian? More specifically, how would you > write an appropriate equation for the heat transfer of a radial system in > Cartesian coordinates? In short, how would you write this: > > > > > > in FiPy? > > > > The CylindricalGrids work fine, except for taking the gradient of the > field, which is an unusual need. > > > > > > > > 3) Ideally, I would move this to a 2D system which is the most > confusing to me. In a 2D system, heat would be transmitted radially while > air would be flowing axially. The air would cool as it passes through the > tube. The diffusion term in the axial direction would be significantly > lower than the convection term. Can I use the same heat transfer equation > you've suggested slightly modified? I would guess the equation to be > something like this: > > > > > > eq = TransientTerm(coeff = rho*Cp,var = T) + > PowerLawConvectionTerm(coeff=rho*Cp*v,var = T) == DiffusionTerm(coeff=k,var > = T) > > > > > > s.t. v = velocity of the air. > > > > Looks OK > > > > > I would also think that I would have to set the value of v to 0 at the > wall and beyond, which means v would be a CellVariable like rho and k. > > > > Actually, v should be a rank-1 FaceVariable. > > > > Anisotropic diffusion is supported in FiPy (see > https://www.ctcms.nist.gov/fipy/examples/diffusion/generated/examples.diffusion.anisotropy.html), > but are you really sure
Re: Question regarding Boundary Condition Implementation
Daniel - Two things going on here: - Your constraint is being applied at the tip of the wedge of the cylindrical coordinates. Since this "face" has no area, that Dirichlet condition doesn't generate any flux in or out of the domain. Offsetting the cylinder slightly from the origin resolves this: mesh = CylindricalGrid1D(Lr = L,nr=nx) + [[L/nx]] I recommend giving some thought as to whether it even makes sense to apply a Dirichlet condition along the centerline of a cylinder. - Because you've introduced coefficients which change with temperature, your equation is now nonlinear and so sweeping is required (2 is apparently enough). At each time step, execute: for sweep in range(2): res = eq.sweep(dt=dt) rho.setValue(rho_l,where = PCM & (T>350) ) k.setValue(k_l,where = PCM & (T>350)) Cp.setValue(Cp_l,where = PCM & (T>350)) With these changes, the evolution seems generally similar to your Grid1D case. - Jon > On Apr 15, 2019, at 1:28 PM, Daniel DeSantis wrote: > > Jonathan, > > I have just a few follow up questions. Namely, regarding to working within > cylindrical coordinates. If I switch the mesh to cylindrical, it seems the > temperature profile never changes, unless I go to extremely small time steps. > That strikes me as odd. What do you think? > > Also, the shifting profile seems to reverse, with the left boundary condition > dropping to ~300K over time while in Cartesian coordinates, the temperature > rises to 600K over time (a result I would expect). I've included a copy of > the work with the changes you've suggested previously. If you look at lines > 23/24, you can toggle the cylindrical grid on and off, and in lines 30-34, > I've included options for the time steps that work for each respective > coordinate system. Do you have any suggestions? > > Finally, in reference to the velocity I would like to incorporate, I > understand the difference between Cell Variables and Face Variables, but why > is the rank of the velocity -1? What does the Face Variable rank indicate? > > Thank you, > Dan > > On Thu, Apr 11, 2019 at 8:14 AM Guyer, Jonathan E. Dr. (Fed) via fipy > wrote: > > > > On Apr 10, 2019, at 4:42 PM, Daniel DeSantis wrote: > > > > 1) Set a different Cp, rho, and k in just the PCM domain when the > > temperature in that domain crosses over a certain melt temperature (350K). > > I tried an if statement inside the solving loop and I think that has > > worked, based on some print testing I did in the loop. I just wanted to get > > your opinion and see if that was the right way to go. > > Seems reasonable, as long as you're calling `.setValue()` and not redefining > Cp, rho, and k. > > In fact, I wouldn't use an `if`: > > >>> for step in range(steps): > ... T.updateOld() > ... eq.solve(dt=dt) > ... rho.setValue(rho_melt, where=PCM & (T > 350)) > > > > 2) I'd like to move this to cylindrical coordinates with the lengths being > > radii. I'm not sure if the cylindrical grid is working though. And if it > > isn't, is there a way to convert the diffusion equation into something that > > would solve the problem in Cartesian? More specifically, how would you > > write an appropriate equation for the heat transfer of a radial system in > > Cartesian coordinates? In short, how would you write this: > > > > in FiPy? > > The CylindricalGrids work fine, except for taking the gradient of the field, > which is an unusual need. > > > > > 3) Ideally, I would move this to a 2D system which is the most confusing to > > me. In a 2D system, heat would be transmitted radially while air would be > > flowing axially. The air would cool as it passes through the tube. The > > diffusion term in the axial direction would be significantly lower than the > > convection term. Can I use the same heat transfer equation you've suggested > > slightly modified? I would guess the equation to be something like this: > > > > eq = TransientTerm(coeff = rho*Cp,var = T) + > > PowerLawConvectionTerm(coeff=rho*Cp*v,var = T) == DiffusionTerm(coeff=k,var > > = T) > > > > s.t. v = velocity of the air. > > Looks OK > > > I would also think that I would have to set the value of v to 0 at the wall > > and beyond, which means v would be a CellVariable like rho and k. > > Actually, v should be a rank-1 FaceVariable. > > Anisotropic diffusion is supported in FiPy (see > https://www.ctcms.nist.gov/fipy/examples/diffusion/generated/examples.diffusion.anisotropy.html), > but are you really sure you have lower axial diffusivity? It seems to me > more likely that transport is simply dominated by axial flow. > > > I also guess this equation would be subject to any changes that would be > > made in question 2, regarding writing a cylindrical equation in a cartesian > > system. What do you think? > > Use a CylindricalGrid and don't change the equation at all. > > > Again, thank you so much for your help and I hope I haven't
Re: Question regarding Boundary Condition Implementation
Jonathan, I have just a few follow up questions. Namely, regarding to working within cylindrical coordinates. If I switch the mesh to cylindrical, it seems the temperature profile never changes, unless I go to extremely small time steps. That strikes me as odd. What do you think? Also, the shifting profile seems to reverse, with the left boundary condition dropping to ~300K over time while in Cartesian coordinates, the temperature rises to 600K over time (a result I would expect). I've included a copy of the work with the changes you've suggested previously. If you look at lines 23/24, you can toggle the cylindrical grid on and off, and in lines 30-34, I've included options for the time steps that work for each respective coordinate system. Do you have any suggestions? Finally, in reference to the velocity I would like to incorporate, I understand the difference between Cell Variables and Face Variables, but why is the rank of the velocity -1? What does the Face Variable rank indicate? Thank you, Dan On Thu, Apr 11, 2019 at 8:14 AM Guyer, Jonathan E. Dr. (Fed) via fipy < fipy@nist.gov> wrote: > > > > On Apr 10, 2019, at 4:42 PM, Daniel DeSantis wrote: > > > > 1) Set a different Cp, rho, and k in just the PCM domain when the > temperature in that domain crosses over a certain melt temperature (350K). > I tried an if statement inside the solving loop and I think that has > worked, based on some print testing I did in the loop. I just wanted to get > your opinion and see if that was the right way to go. > > Seems reasonable, as long as you're calling `.setValue()` and not > redefining Cp, rho, and k. > > In fact, I wouldn't use an `if`: > > >>> for step in range(steps): > ... T.updateOld() > ... eq.solve(dt=dt) > ... rho.setValue(rho_melt, where=PCM & (T > 350)) > > > > 2) I'd like to move this to cylindrical coordinates with the lengths > being radii. I'm not sure if the cylindrical grid is working though. And if > it isn't, is there a way to convert the diffusion equation into something > that would solve the problem in Cartesian? More specifically, how would you > write an appropriate equation for the heat transfer of a radial system in > Cartesian coordinates? In short, how would you write this: > > > > in FiPy? > > The CylindricalGrids work fine, except for taking the gradient of the > field, which is an unusual need. > > > > > 3) Ideally, I would move this to a 2D system which is the most confusing > to me. In a 2D system, heat would be transmitted radially while air would > be flowing axially. The air would cool as it passes through the tube. The > diffusion term in the axial direction would be significantly lower than the > convection term. Can I use the same heat transfer equation you've suggested > slightly modified? I would guess the equation to be something like this: > > > > eq = TransientTerm(coeff = rho*Cp,var = T) + > PowerLawConvectionTerm(coeff=rho*Cp*v,var = T) == DiffusionTerm(coeff=k,var > = T) > > > > s.t. v = velocity of the air. > > Looks OK > > > I would also think that I would have to set the value of v to 0 at the > wall and beyond, which means v would be a CellVariable like rho and k. > > Actually, v should be a rank-1 FaceVariable. > > Anisotropic diffusion is supported in FiPy (see > https://www.ctcms.nist.gov/fipy/examples/diffusion/generated/examples.diffusion.anisotropy.html), > but are you really sure you have lower axial diffusivity? It seems to me > more likely that transport is simply dominated by axial flow. > > > I also guess this equation would be subject to any changes that would be > made in question 2, regarding writing a cylindrical equation in a cartesian > system. What do you think? > > Use a CylindricalGrid and don't change the equation at all. > > > Again, thank you so much for your help and I hope I haven't taken up too > much of your time. > > Happy to help. > > > ___ > fipy mailing list > fipy@nist.gov > http://www.ctcms.nist.gov/fipy > [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ] > -- Daniel DeSantis # -*- coding: utf-8 -*- """ Created on Tue Apr 2 15:14:40 2019 @author: ddesantis """ from fipy import * import numpy as np import matplotlib.pyplot as plt #%% - Mesh nx = 100 L = 0.061 #m Ly = 1 #m L1 = 0.035 # m, Length of air channel t_w = 1./1000. # m, wall thickness L2 = L1+t_w # m, interface between PCM and wall L3 = L # m, Total Length mesh = Grid1D(Lx=L,nx=nx) #mesh = CylindricalGrid1D(Lr = L,nr=nx) x = mesh.cellCenters[0] #%% - Steps # Cartesian coordinates time step dt = 0.25 #s,time step for sweeps # Cylindrical coordinates time step #dt = 0.5 #s,time step for sweeps steps = 150 #%% Universal Constants MW = 2.02#g H2 mol^-1 R = 8.314e-3# kj mol^-1 K^-1 #%% - air Properties T_in = 600. # K P_in = 1. # bar
Re: Question regarding Boundary Condition Implementation
> On Apr 10, 2019, at 4:42 PM, Daniel DeSantis wrote: > > 1) Set a different Cp, rho, and k in just the PCM domain when the temperature > in that domain crosses over a certain melt temperature (350K). I tried an if > statement inside the solving loop and I think that has worked, based on some > print testing I did in the loop. I just wanted to get your opinion and see if > that was the right way to go. Seems reasonable, as long as you're calling `.setValue()` and not redefining Cp, rho, and k. In fact, I wouldn't use an `if`: >>> for step in range(steps): ... T.updateOld() ... eq.solve(dt=dt) ... rho.setValue(rho_melt, where=PCM & (T > 350)) > 2) I'd like to move this to cylindrical coordinates with the lengths being > radii. I'm not sure if the cylindrical grid is working though. And if it > isn't, is there a way to convert the diffusion equation into something that > would solve the problem in Cartesian? More specifically, how would you write > an appropriate equation for the heat transfer of a radial system in Cartesian > coordinates? In short, how would you write this: > > in FiPy? The CylindricalGrids work fine, except for taking the gradient of the field, which is an unusual need. > > 3) Ideally, I would move this to a 2D system which is the most confusing to > me. In a 2D system, heat would be transmitted radially while air would be > flowing axially. The air would cool as it passes through the tube. The > diffusion term in the axial direction would be significantly lower than the > convection term. Can I use the same heat transfer equation you've suggested > slightly modified? I would guess the equation to be something like this: > > eq = TransientTerm(coeff = rho*Cp,var = T) + > PowerLawConvectionTerm(coeff=rho*Cp*v,var = T) == DiffusionTerm(coeff=k,var = > T) > > s.t. v = velocity of the air. Looks OK > I would also think that I would have to set the value of v to 0 at the wall > and beyond, which means v would be a CellVariable like rho and k. Actually, v should be a rank-1 FaceVariable. Anisotropic diffusion is supported in FiPy (see https://www.ctcms.nist.gov/fipy/examples/diffusion/generated/examples.diffusion.anisotropy.html), but are you really sure you have lower axial diffusivity? It seems to me more likely that transport is simply dominated by axial flow. > I also guess this equation would be subject to any changes that would be made > in question 2, regarding writing a cylindrical equation in a cartesian > system. What do you think? Use a CylindricalGrid and don't change the equation at all. > Again, thank you so much for your help and I hope I haven't taken up too much > of your time. Happy to help. ___ fipy mailing list fipy@nist.gov http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
Re: Question regarding Boundary Condition Implementation
Jonathan, Thank you very much for you help. It wouldn't be the first time I've over complicated a system. I understand how you set up the system a little bit better now. I had set things up the way I had before because I had expected to expand the system a bit. I was hoping to do the following: 1) Set a different Cp, rho, and k in just the PCM domain when the temperature in that domain crosses over a certain melt temperature (350K). I tried an if statement inside the solving loop and I think that has worked, based on some print testing I did in the loop. I just wanted to get your opinion and see if that was the right way to go. 2) I'd like to move this to cylindrical coordinates with the lengths being radii. I'm not sure if the cylindrical grid is working though. And if it isn't, is there a way to convert the diffusion equation into something that would solve the problem in Cartesian? More specifically, how would you write an appropriate equation for the heat transfer of a radial system in Cartesian coordinates? In short, how would you write this: [image: image.png] in FiPy? 3) Ideally, I would move this to a 2D system which is the most confusing to me. In a 2D system, heat would be transmitted radially while air would be flowing axially. The air would cool as it passes through the tube. The diffusion term in the axial direction would be significantly lower than the convection term. Can I use the same heat transfer equation you've suggested slightly modified? I would guess the equation to be something like this: eq = TransientTerm(coeff = rho*Cp,var = T) + PowerLawConvectionTerm(coeff=rho*Cp*v,var = T) == DiffusionTerm(coeff=k,var = T) s.t. v = velocity of the air. I would also think that I would have to set the value of v to 0 at the wall and beyond, which means v would be a CellVariable like rho and k. I also guess this equation would be subject to any changes that would be made in question 2, regarding writing a cylindrical equation in a cartesian system. What do you think? Again, thank you so much for your help and I hope I haven't taken up too much of your time. Thank you, Dan On Mon, Apr 8, 2019 at 5:37 PM Guyer, Jonathan E. Dr. (Fed) via fipy < fipy@nist.gov> wrote: > You're over-complicating things trying to have three different > temperatures. There's one temperature in three different domains, with > different material properties in each domain. FiPy (and physics) takes care > of matching fluxes for you at the internal boundaries. > > > The changes I made are at > https://gist.github.com/guyer/d86ee6f085e9832df781286099a73800/revisions > > The primary changes I made: > > - I defined three domains, `air`, `wall`, and `PCM`, and then defined `T`, > `rho`, `Cp`, and `k` to have different values in those different domains. > > - sweeps are for getting better convergence of non-linear equations. > Your equations are linear and from your usage you want time steps, not > sweeps. > See: > https://www.ctcms.nist.gov/fipy/documentation/FAQ.html#iterations-timesteps-and-sweeps-oh-my > > - Internal boundary conditions need to be defined very differently (see: > https://www.ctcms.nist.gov/fipy/documentation/USAGE.html#applying-internal-boundary-conditions > ) > Importantly, you don't need internal boundary conditions. > > - All field variables and coefficients should be floats, not integers. > > > Note: Your wall thickness `tw = 0.001` is smaller than your grid spacing > `dx = 0.061 / 50 = 0.00122`, your your numerics will be terrible. Increase > your wall thickness or decrease your grid size. > > > On Apr 8, 2019, at 2:44 PM, Daniel DeSantis wrote: > > > > Hello, > > > > I am trying to model heat diffusion through a wall using FiPy. > Essentially, hot air on one side of the wall would heat the wall, and the > wall would subsequently heat a phase change material on the other side (See > picture below). The eventual goal would be to have the phase change > material properties change as the heat increases. I am having a few issues > with setting boundary conditions however. > > > > The goal of the model would be to have the heat flux at the wall be > equal to the heat flux of the fluid on each respective side. Further, the > heat flux of one of the fluids (a phase change material) would be 0 at the > side not in contact with the wall. (I've included a copy of my FiPy code as > is, and a series of equations and boundary conditions below). > > > > Right now, the temperature for the air and wall are staying constant, > while the phase change material temperature drops. I think how I've defined > the boundary conditions are to blame, but I'm not sure what I'm doing > wrong. Any help here would be greatly appreciated. > > > > Thank you, > > > > -- > > Daniel DeSantis > > > > > > > > > > > > > > > ___ > > fipy mailing list > > fipy@nist.gov > > http://www.ctcms.nist.gov/fipy > > [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ] > > >
Re: Question regarding Boundary Condition Implementation
You're over-complicating things trying to have three different temperatures. There's one temperature in three different domains, with different material properties in each domain. FiPy (and physics) takes care of matching fluxes for you at the internal boundaries. The changes I made are at https://gist.github.com/guyer/d86ee6f085e9832df781286099a73800/revisions The primary changes I made: - I defined three domains, `air`, `wall`, and `PCM`, and then defined `T`, `rho`, `Cp`, and `k` to have different values in those different domains. - sweeps are for getting better convergence of non-linear equations. Your equations are linear and from your usage you want time steps, not sweeps. See: https://www.ctcms.nist.gov/fipy/documentation/FAQ.html#iterations-timesteps-and-sweeps-oh-my - Internal boundary conditions need to be defined very differently (see: https://www.ctcms.nist.gov/fipy/documentation/USAGE.html#applying-internal-boundary-conditions) Importantly, you don't need internal boundary conditions. - All field variables and coefficients should be floats, not integers. Note: Your wall thickness `tw = 0.001` is smaller than your grid spacing `dx = 0.061 / 50 = 0.00122`, your your numerics will be terrible. Increase your wall thickness or decrease your grid size. > On Apr 8, 2019, at 2:44 PM, Daniel DeSantis wrote: > > Hello, > > I am trying to model heat diffusion through a wall using FiPy. Essentially, > hot air on one side of the wall would heat the wall, and the wall would > subsequently heat a phase change material on the other side (See picture > below). The eventual goal would be to have the phase change material > properties change as the heat increases. I am having a few issues with > setting boundary conditions however. > > The goal of the model would be to have the heat flux at the wall be equal to > the heat flux of the fluid on each respective side. Further, the heat flux of > one of the fluids (a phase change material) would be 0 at the side not in > contact with the wall. (I've included a copy of my FiPy code as is, and a > series of equations and boundary conditions below). > > Right now, the temperature for the air and wall are staying constant, while > the phase change material temperature drops. I think how I've defined the > boundary conditions are to blame, but I'm not sure what I'm doing wrong. Any > help here would be greatly appreciated. > > Thank you, > > -- > Daniel DeSantis > > > > > > > ___ > 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 ]
Question regarding Boundary Condition Implementation
Hello, I am trying to model heat diffusion through a wall using FiPy. Essentially, hot air on one side of the wall would heat the wall, and the wall would subsequently heat a phase change material on the other side (See picture below). The eventual goal would be to have the phase change material properties change as the heat increases. I am having a few issues with setting boundary conditions however. The goal of the model would be to have the heat flux at the wall be equal to the heat flux of the fluid on each respective side. Further, the heat flux of one of the fluids (a phase change material) would be 0 at the side not in contact with the wall. (I've included a copy of my FiPy code as is, and a series of equations and boundary conditions below). Right now, the temperature for the air and wall are staying constant, while the phase change material temperature drops. I think how I've defined the boundary conditions are to blame, but I'm not sure what I'm doing wrong. Any help here would be greatly appreciated. Thank you, -- Daniel DeSantis [image: image.png] [image: image.png] [image: image.png] [image: image.png] PCM_thermal_v0-SendOut.py Description: Binary data ___ fipy mailing list fipy@nist.gov http://www.ctcms.nist.gov/fipy [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
