Thanks Daniel and Jonathan,
All of this information helps very much.
It is now working properly, so here is the code.
In case it is useful to anyone.
It is basic, but I figure its nice to share.
cellSize=0.05
radius=1.
from fipy import *
mesh = Gmsh2D('''cellSize = %(cellSize)g;
radius = %(radius)g;
Point(1) = {0, 0, 0, cellSize};
Point(2) = {radius, 0, 0, cellSize};
Point(3) = {0, radius, 0, cellSize};
Point(5) = {-radius, 0, 0, cellSize};
Point(6) = {0, -radius, 0, cellSize};
Circle(1) = {2, 1, 3};
Circle(2) = {3, 1, 5};
Circle(3) = {5, 1, 6};
Circle(4) = {6, 1, 2};
Line Loop(1) = {1, 2, 3, 4};
Point(200) = {radius/3, 0, 0, cellSize};
Point(300) = {0, radius/3, 0, cellSize};
Point(500) = {-radius/3, 0, 0, cellSize};
Point(600) = {0, -radius/3, 0, cellSize};
Circle(100) = {200, 1, 300};
Circle(200) = {300, 1, 500};
Circle(300) = {500, 1, 600};
Circle(400) = {600, 1, 200};
Line Loop(100) = {100, 200, 300, 400};
Plane Surface(1)={1,100};
''' % locals())
phi=CellVariable(name="solution",mesh=mesh,value=0.)
D=1
eq=TransientTerm()==DiffusionTerm(coeff=D)
X,Y=mesh.faceCenters
phi.constrain(1,mesh.exteriorFaces)
timeStepDuration=10*0.9*cellSize**2/(2*D)
steps=10
for step in range(steps):
eq.solve(var=phi,dt=timeStepDuration)
TSVViewer(vars=(phi)).plot(filename="diffusion.annulus."+str(step)+".tsv”)
Thanks!
Kyle.
On Nov 19, 2014, at 2:16 PM, Guyer, Jonathan E. Dr. <[email protected]>
wrote:
> Just to clarify: as far as FiPy is concerned, "interiorFaces" have a cell on
> either side of them, whereas "exteriorFaces" have a cell on only one side;
> the other side is "outside" the mesh. It doesn't matter whether the mesh has
> simple topological connectivity (square, sphere, etc.) or complex topological
> connectivity (torus). "interiorFaces" separate "mesh" from "more mesh";
> "exteriorFaces" separate "mesh" from "not mesh".
>
>
> On Nov 19, 2014, at 10:53 AM, Daniel Wheeler <[email protected]>
> wrote:
>
>> Hi Kyle,
>>
>> #Boundary Conditions
>> phi.constrain(X,mesh.exteriorFaces)
>> phi.constrain(X,mesh.interiorFaces)
>>
>> The above is the problem. You are constraining the internal faces, which
>> makes no sense in FiPy. I am not even sure how FiPy behaves when that
>> constraint is added. However, I assume that is not what you want to do. If
>> you remove that constraint the result seems to look nice. BTW,
>> "exteriorFaces" refers to both the inner and outer external faces of the
>> annulus. Perhaps you thought "interiorFaces" referred to the interior
>> annulus faces that are external to the mesh.
>>
>> Also, you don't need to import the results into Mathematica. You can just
>> use "viewer = fipy.Viewer(phi); viewer.plot()" to look at the results.
>>
>> Cheers,
>>
>> Daniel
>>
>> --
>> Daniel Wheeler
>> _______________________________________________
>> fipy mailing list
>> [email protected]
>> http://www.ctcms.nist.gov/fipy
>> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
>
>
> _______________________________________________
> fipy mailing list
> [email protected]
> http://www.ctcms.nist.gov/fipy
> [ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
_______________________________________________
fipy mailing list
[email protected]
http://www.ctcms.nist.gov/fipy
[ NIST internal ONLY: https://email.nist.gov/mailman/listinfo/fipy ]
